CN114170524A

CN114170524A - Single-day tidal flat digital terrain construction method based on high-time-space satellite data fusion

Info

Publication number: CN114170524A
Application number: CN202111262990.2A
Authority: CN
Inventors: 高文丽; 沈芳; 谭凯; 张卫国; 刘权兴; 葛建忠
Original assignee: East China Normal University
Current assignee: East China Normal University
Priority date: 2021-10-28
Filing date: 2021-10-28
Publication date: 2022-03-11

Abstract

The invention discloses a single-day tidal flat digital terrain construction method for high time-space satellite data fusion, which fuses GF-1/WFV satellite data and geostationary satellite GOCI data by using a STARFM time-space data fusion method to generate hourly images with the spatial resolution of 16 m; the tidal water level at the image transit time is obtained by combining the hydrodynamic model, the daily terrain elevation of the tidal flat is obtained, and the remote sensing dynamic monitoring of the tidal flat terrain is realized. Through the comprehensive application of high-resolution satellite observation and dynamic simulation, the difficulty of dynamic monitoring of tidal flat terrain is solved, and a new method for dynamic monitoring of silt tidal flat terrain is established.

Description

Single-day tidal flat digital terrain construction method based on high-time-space satellite data fusion

Technical Field

The invention belongs to the technical field of remote sensing construction of tidal flat terrain Digital Elevation models (Digital Elevation models DEMs), and particularly relates to a remote sensing method for generating single-day tidal flat terrain Digital elevations by fusing high-time-high spatial resolution satellite data.

Background

The tidal flat is a shoal with exposed intertidal zones, contains precious wetland biological resources and land resources, is influenced by human engineering, sea level rise and the like, and silts and erosion of the shoal are converted along with the tidal flat, so that the terrain elevation is changed rapidly, and certain threats are caused to the tidal flat wetland resources and the land resources. Therefore, the method has important significance for timely and real-timely grasping the rapid change of the tidal flat terrain elevation. In general, conventional topographic surveying techniques (e.g., leveling, total station, etc.) can provide very accurate tidal flat topographic elevation measurements, but the survey space coverage and work efficiency are very low, and a lot of manpower and financial resources are consumed, and furthermore, tidal flat wetland, mudflat, mud and marsh environments are very difficult to measure under the influence of tidal water fluctuation. The satellite remote sensing technology has the advantages of space-time continuous and large-area synchronous observation, short measurement period, low cost, non-contact type and the like, and realizes the tidal flat terrain measurement by constructing the DEM through the instantaneous water line and the synchronous water level information of the satellite. Generally, optical remote sensing images are richer in data and easier to obtain compared with radar images and stereopair images. Therefore, it is currently most widespread and practical to construct tidal flat digital terrain elevation by using optical remote sensing image data.

Due to tidal power factors, coastal zone tidal beaches are in dynamic change at any moment, but for the remote sensing technology, the spatial distribution range of the coastal zone silt tidal beaches is not large, and a satellite image with high resolution is required. The number of time-frequency low high-resolution satellite images is extremely limited, and in addition, due to the influence of cloud covers in coastal zone areas, the DEM is difficult to construct under the condition that a time sequence series image data set with high-low tide level changes is met in a short time, most scholars adopt different data sources to achieve the condition, the shortest time is shortened to within 1 year, within 0.5 year and the like, but the rapid dynamic change of the wetland in the intertidal zone cannot be really solved.

Disclosure of Invention

The invention aims to provide a method for constructing a single-day tidal flat terrain DEM (digital elevation model) by data space-time fusion of hourly revisit frequency and low spatial resolution (hundred meter level) remote sensing images of a static orbit satellite and multi-day revisit frequency and high spatial resolution (ten meter level) remote sensing images of an polar orbit satellite. The method has the advantages that a single-day tidal flat terrain digital elevation model is constructed, dynamic and rapid change monitoring of tidal flat terrain elevation is achieved, and the major requirements of real-time dynamic monitoring of tidal flat wetland resources and sustainable development of land resources in national economy and ecological civilization construction are met.

The specific technical scheme for realizing the purpose of the invention is as follows:

a single-day tidal flat digital terrain construction method for high-time-space satellite data fusion is characterized by comprising the following steps: the method comprises the following specific steps:

step 1, fusing high-space-time resolution images, and generating hourly images with the spatial resolution of 16m by adopting a space-time adaptive reflection fusion model STARFM;

step 2, acquiring a tide level by using a tide power model TPXO8, and acquiring tide data by using a tide power model TPXO 8;

step 3, based on the water line extraction of the BP neural network, accurately classifying the muddy tidal wetland and the water body by adopting a BP neural network supervision classification method, and further extracting the water line;

and 4, constructing the DTM, and constructing the tidal flat terrain elevation by combining the tidal level and the water line data.

The step 1 specifically comprises:

step 1.1, firstly, introducing information of adjacent and similar pixels in a window, firstly resampling a GOCI image after geometric correction and geographic registration to a GF-1/WFV image with 16m resolution by adopting a STARFM algorithm, and setting pixels (x)_i,y_j) Type of ground coverage and system error over time t₀And image date t_kIf not, then:

GF(x_i,y_j,t₀)＝G(x_i,y_j,t₀)+GF(x_i,y_j,t_k)-G(x_i,y_j,t_k) (1)

in the formula: GF and G represent GF-1/WFV and GOCI pictures, respectively, (x)_i,y_j) Is the pixel position, t₀To predict time, t_kTime of an existing image pair;

step 1.2, predicting the central pixel value by utilizing a sliding window and a weight function

The STARFM algorithm assumes that the errors of the GF-1/WFV and GOCI imagery are constant for the same ground target, and that the errors are systematic in a short time, relying only on the characteristics of the individual pixels; therefore, the error of each pixel is estimated by using a pair of GF-1/WFV-GOCI images at a certain time, and then the high spatial resolution GF-1/WFV images are predicted by using the GOCI images at other times;

the STARFM algorithm first resamples the geometrically corrected and geographically registered GOCI image to a GF-1/WFV image of 16m resolution, sets the pixels (x)_i,y_j) Type of ground coverage and system error over time t₀And image date t_kIf not, then:

GF(x_i,y_j,t₀)＝G(x_i,y_j,t₀)+GF(x_i,y_j,t_k)-G(x_i,y_j,t_k) (1)

however, the GOCI image is usually a mixed pixel under the GF-1/WFV 16m resolution, or within the forecast date, the tidal bank water side line is often changed due to the influence of tide; the predicted time t is calculated using the weighting function in the STARFM algorithm by introducing additional information of neighboring pixels₀Surface reflectance of center pixel:

in the formula: ω is the size of the search window, (x)_ω/2,y_ω/2) Is the center pixel position of the search window; weight W_ijkDetermining a contribution of each neighboring pixel to the estimated reflectivity of the central pixel; w_ijkFrom spectral differences, time differences, central pixel position (x)_ω/2,y_ω/2) And candidate pixel location (x)_i,y_j) The distance between them is determined;

the spectral difference between GF-1/WFV and GOCI data at a given location is:

S_ijk＝|GF(x_i,y_j,t_k)-G(x_i,y_j,t_k)| (3)

the difference between the GF-1/WFV observation value and the average value of the adjacent pixels under low spatial resolution is approximately obtained by formula (3); s_ijkA smaller value of (a) indicates that the pixel has a spectral characteristic closer to the average of its surrounding pixels, i.e., the variation of the high resolution pixel is comparable to the average variation of the surrounding pixels, where the reflectance of the pixel is assigned a higher weight in equation (2); the time difference between the input GOCI data and the predicted GOCI data is:

T_ijk＝|G(x_i,y_j,t_k)-G(x_i,y_j,t₀)| (4)

T_ijkthe variation between prediction and acquisition times is measured, the smaller T_ijkMeaning the time t_kAnd t₀The change between is smaller and the pixel should be assigned a higher weight; at the date t_kTemporal center pixel position (x)_ω/2,y_ω/2) And candidate pixel location (x)_i,y_j) The distance between them is:

d_ijkthe spatial distance between the center predicted pixel and the surrounding spectrally similar candidate pixel is assigned a higher weight for its candidate pixel for closer pixels.

The step 2 specifically comprises:

step 2.1, the tidal height of the tidal model TPXO8 is given in the form of complex amplitude, so the partial tide at time t for a single component of frequency f at location x is given by the following equation:

h(t,x)＝pu(t,x)·Re[h(x)exp{i[f(t-t₀)+V₀(t₀)+ph(t,x)]}] (6)

wherein V₀(t₀) Is t₀The astronomical parameters at time, pu (t, x) and ph (t, x) are node correction quantities; then the amplitude is | h | and the phase is tan^-1(-IM (h)/RE (h)), where IM and RE represent imaginary Number part and Real Number part of complex Number, i x i ═ 1;

and 2.2, obtaining Matlab toolbox codes from the OSU-TPXO Tide Models (https:// www.tpxo.net/home), and inputting the accurate position coordinates and time information into the TPXO8 model to obtain Tide data.

The step 3 specifically includes:

step 3.1, selecting two types of interested areas of the intertidal zone wetland and the water body as training samples, and taking the training samples as an input layer for supervision and classification;

step 3.2, utilizing a BP neural network based on ENVI5.5 (trial version) to supervise and classify the images;

and 3.3, outputting the supervised classification result and the separation degree of the sample, wherein the separation degree of the sample is larger if the difference between the two classification areas is obvious, the separation degree is better, the separation degree is maximum 2.0, and a parting line between the two classification areas is a water line.

The step 4 specifically includes:

step 4.1, constructing a Delaunay (triangulation method) irregular triangulation network, and performing subdivision by using the Delaunay irregular triangulation;

and 4.2, converting TIN (triangulation) into grids, and interpolating the numerical values into a grid DTM (digital terrain elevation) to obtain the digital terrain elevation.

The invention has the beneficial effects that:

compared with the prior art, the method has the advantages that the time interval of the images is shortened to 1h, the spatial resolution is 16m, and the monitoring of the DTM of the tidal flat is facilitated. The tidal flat wetland resource real-time dynamic monitoring is realized, and a foundation and a premise are provided for the protection and treatment of the estuary coastal zone.

Drawings

FIG. 1 is a flow chart of the present invention;

FIG. 2 is a STARFM spatiotemporal fusion flow diagram;

fig. 3 is a TPXO8 model interface diagram;

FIG. 4 is a flow chart of the waterside line extraction;

FIG. 5 is a flow chart of a DTM construction using a water sideline;

FIG. 6 is a DTM verification graph constructed from Lianxing harbor measured tide level and TPXO8 model tide level;

FIG. 7 is a Chongming east beach DTM and difference plot across a one-year scale;

FIG. 8 is a comparison of water lines at 500m and 16m spatial resolution;

FIG. 9 is a DTM_LIDAnd annual DTMON-fusion (DTM)_LID-DTMnon-fusion).

Detailed Description

Referring to fig. 1, the present invention performs data fusion of GF-1/WFV (spatial resolution: 16 m; revisit frequency: 4 days) and GOCI (spatial resolution: 500 m; revisit frequency: 1 hour) using a spatio-temporal adaptive reflection fusion model (STARFM), and generates an hourly video with a spatial resolution of 16 m. Extracting a water line by using a BP (Back propagation) neural network supervision and classification method, and constructing a digital ground model (DTM) of the tidal flat of the Yangtze river estuary by combining with the tide level data acquired by the TPXO8 tide model; the method specifically comprises the following steps:

step 1: high spatial-temporal resolution image fusion

In order to obtain high-precision intertidal wetland terrain elevation by using a water line extraction method, a satellite image is required to have higher spatial resolution so as to meet the precision requirement of DTM. Meanwhile, for the water line extraction, the revisit period of the satellite images is required to be short enough to capture the rapid change dynamics of the intertidal zone wetland elevation. Aiming at the problem of insufficient space-time resolution of remote sensing images, the STARFM fusion method is used for fusing GF-1 data with high spatial resolution and GOCI data. The model considers the difference of space and the difference of time in the fusion process, and is one of the most widely applied space-time fusion models at present. Compared with other fusion models, STARFM is more suitable for regions dominated by spatio-temporal changes, and is very suitable for the study of dynamically changing tidal flat DTM. Therefore, STARFM is adopted to fuse GF-1 and GOCI data, and the space-time resolution of the remote sensing image is further improved.

STARFM uses a sliding window and a weight function to predict the central pixel value of the window by introducing information of neighboring pixels and similar pixels in the window. The algorithm assumes that the errors of the GF-1 and GOCI shots are invariant when observing the same ground target. The error depends on the characteristics of the individual pixels and is systematic over a short time interval. Therefore, if one image pair of GF-1 and GOCI can be obtained as a reference at a certain time, an error of each pixel in the image can be estimated. And then predicting the fine resolution image by using the GOCI images at other time instants.

Step 1.1: the STARFM algorithm first resamples the geometry-corrected and geo-registered GOCI image onto a GF-1 image of 16m resolution, assuming pixels (x)_i,y_j) Type of ground coverage and system error over time t₀And image date t_kIf not, then:

GF(x_i,y_j,t₀)＝G(x_i,y_j,t₀)+GF(x_i,y_j,t_k)-G(x_i,y_j,t_k) (1)

in the formula: GF and G represent GF-1/WFV and GOCI pictures, respectively, (x)_i,y_j) Is the pixel position, t₀To predict time, t_kThe existing image is compared to time.

Step 1.2: because GOCI images are usually mixed pixels at GF-116 m resolution, or because of tides and the like within the forecast date, tidal shoal water lines are often changed. Thus STARFM uses a sliding window and a weighting function to predict the value of the center pel within the sliding window by introducing information of neighboring and similar pels within the window. By introducing additional information on neighboring pixels, the STARFM algorithm calculates the surface reflectance of the pixel in the center of the predicted date using a weighting function:

in the formula: ω is the size of the search window, (x)_ω/2,y_ω/2) Is the center pixel position of the search window. Weight W_ijkThe contribution of each neighboring pixel to the estimated reflectivity of the central pixel is determined. W_ijkFrom spectral differences, time differences, central pixel position (x)_ω/2,y_ω/2) And candidate pixel location (x)_i,y_j) The distance between them, etc.

The spectral difference between GF-1 and GOCI data at a given location is:

S_ijk＝|GF(x_i,y_j,t_k)-G(x_i,y_j,t_k)| (3)

equation (3) approximately measures the difference between the GF-1 observation and the average of the neighboring pixels at coarse resolution. S_ijkA smaller value of (a) indicates that the pixel has a spectral characteristic closer to the average of the surrounding pixels, i.e. the variation of the high resolution pixel should be comparable to the average variation of the surrounding pixels. Therefore, the reflectance of the pixel should be assigned a higher weight in equation (2). Input GOCI dataAnd the predicted GOCI data is:

T_ijk＝|G(x_i,y_j,t_k)-G(x_i,y_j,t₀)| (4)

T_ijkthe variation between prediction and acquisition time is measured. Smaller T_ijkMeaning times t and t₀The change between is smaller and the pixel should be assigned a higher weight. At the date t_kTemporal center pixel position (x)_ω/2,y_ω/2) And candidate pixel location (x)_i,y_j) The distance between them is:

d_ijkthe spatial distance between the center predicted pixel and the surrounding spectrally similar candidate pixel is predicted. Spatial similarity is generally better for closer pixels, and therefore, closer candidates should be assigned higher weights. A flow chart of the process is shown in fig. 2. Step 2: TPXO8 tidal power model acquisition tidal level

Step 2.1: the tidal height of TPXO8 is given in complex amplitude, so the partial tide at time t for a single component of frequency f at location x is given by the following equation:

h(t,x)＝pu(t,x)·Re[h(x)exp{i[f(t-t₀)+V₀(t₀)+ph(t,x)1}] (6)

wherein V₀(t₀) Is t₀The astronomical parameters at time, pu (t, x) and ph (t, x) are the node correction quantities. Then the amplitude is | h | and the phase is tan^-1(-IM (h)/RE (h)) IM and RE represent the imaginary Number part and the Real Number part of the complex Number i ═ 1, respectively.

Step 2.2: matlab kit code was obtained from OSU-TPXO Tide Models (https:// www.tpxo.net/home). The precise position coordinates and time information are input into the TPXO8 model to obtain tidal level data. The input, simulation and output interfaces are shown in fig. 3.

Step 3, water line extraction and DTM construction based on BP neural network

Considering that a data source is a passive non-stereo optical remote sensing image, a water line detection method is adopted to construct the DTM. In order to automatically extract the water sideline, a supervision and classification method based on a BP neural network in ENVI5.5 (trial edition) is adopted, which consists of hidden layers among an input layer, an output layer and an input-output layer, and respectively corresponds to a selected intertidal wetland and water body supervision sample, a supervision and classification result (intertidal wetland and water body) and a neural network supervision and classification process.

Step 3.1: firstly, two types of interesting areas of the intertidal zone wetland and the water body are selected as training samples.

Step 3.2: then, images were supervised classified using a BP neural network based on ENVI5.5 (trial version). The results show that the error is gradually reduced and the separation degree of the final sample is more than 1.9. Typically, if the difference between the two classification regions is significant, the maximum separability of the sample is 2.0. In the silt tidal flat area, the spectrum difference between the intertidal zone wetland and the water body in the visible light range is not obvious, so that if the separability of the sample is more than 1.9, the differentiation is qualified. The intertidal wetland and the water body are classified, and the boundary between the intertidal wetland and the water body is used as a water boundary line.

And 4, step 4: construction of DTM

Based on arcgis10.1, DTM is constructed with a water line (fig. 5). The DTM is generated by the TIN method of constructing a Delaunay irregular triangulation network.

Step 4.1: firstly, data are generated by utilizing Delaunay irregular triangulation, and the specific operation is realized in ArcGIS 10.1. First, in the 3D analysis extension module of arcgis10.1, TIN is generated by selecting Create TIN.

Step 4.2, interpolating the data into a grid DTM, selecting a TIN to rate in a 3D analysis expansion module of arcgis10.1, converting the TIN into a grid to obtain corresponding DTM data, and an operation flow chart is shown in fig. 5.

Examples

The embodiment realizes dynamic monitoring of tidal flat terrain through a space-time fusion model, a remote sensing image and a tidal model. And performing data fusion on GF-1/WFV and GOCI by using a STARFM space-time adaptive reflection fusion model to generate an hourly image with the spatial resolution of 16 meters. A water line is extracted by using a BP (Back propagation) neural network supervision and classification method, and a Yangtze river mouth tidal flat DTM is constructed by combining tide level data acquired by a TPXO8 tide model.

a. Elevation verification of the Yangtze river mouth Chongming northern beach: with DTM_LiD(contemporaneous unmanned aerial vehicle-mounted LiDAR data) as a ground truth value, and the DTM constructed by the geostationary satellite is verified. Satellite DTM data is based on 2019 annual hourly time series satellite fusion images of 8, 17, 8 and 9 months, in order to verify DTM without actually measured tide level data^msat accuracy using the data of the station of the sea level (DTM) of FIG. 6a^ssat) and the model tide level Data (DTM) of FIG. 6bTPXO8^msat) are respectively assigned as the waterside line elevations to construct satellite DTMs, wherein fig. 6c and 6d are statistical deviation values corresponding to fig. 6a and 6b, and the root mean square errors thereof are 0.22m and 0.16m, respectively. Indicating that the DTM accuracy built using the tide level modeled with TPXO8 is acceptable, it is reasonable to assign a water line using the tide level data of the TPXO8 model in areas without tide stations.

b. Annual change of the landforms of the Yangtze river mouth Chongming east beach: the annual change of the Chongming east beach DTM is researched by the method, 6 time-by-time fusion images of each day of 12-18 days in 2018 and 12-13 days in 2019 are generated by a STARFM space-time fusion method, and two groups of day DTM (shown in FIGS. 7a and 7b) corresponding to 2018 and 2019 are respectively constructed. It was found that there were some areas of greater variation (north in fig. 7 c) in chongming east beach in 2018 and 2019, with significant overall variation, with an average elevation variation of 0.21m and a maximum of 0.98 m. This means that the period of the time-series images should be as short as possible, which is compatible with the method proposed by the present invention.

c. The effect of spatial and temporal resolution on tidal flat terrain elevation: to quantify the effect of spatial resolution on DTM accuracy, high spatial resolution video (GF-1/WFV video, 16m) and low spatial resolution video (GOCI video, 500m) were selected for horizontal line extraction, respectively, for the same region. The maximum horizontal distance between the water line lines of the two data extractions was found to be 372.9m, with an average distance of 134.6 m. If the slope of the tidal flat is assumed to be 0.15 °, the error in estimating the vertical elevation is about 0.35m, and therefore the spatial resolution of the satellite images has a significant impact on the DTM accuracy.

In order to investigate the influence of the time scale of time series image acquisition on the DTM (DTMON-fusion), 6 GF-1/WFV images are used as an experimental area in Chongming North beach, and the DTM (DTMON-fusion) is constructed in a time scale of one year (24 months and 24 days in 2019 to 13 months and 12 months in 2019). The root mean square error is verified to be 0.27m and is larger than 0.16m of daily DTM, which shows that the error of the DTMON-fusion contains the offset caused by the change of the terrain elevation per se within one year.

In fig. 8, (a) is a position diagram, (b) is a water line at a spatial resolution of 500m, and (c) is a water line at a spatial resolution of 16 m. As can be seen from fig. 8(b), the extraction of the water boundary in the GOCI image is affected by the mixed pixels due to the low resolution, resulting in a large deviation of the position of the water boundary, compared to the water boundary extracted in the GF-1/WFV image which is less affected. The maximum horizontal distance between the two waterlines is 372.9m, and the average distance between the two straight lines is 134.6m (the area between the two straight lines divided by the straight line distance between the two ends of the straight line) calculated by ArcGISI 10.1. In the local area, the gradient of the intertidal zone wetland is 0.15 degrees, the error in the horizontal distance is 134.6m, the error in the vertical elevation can be calculated to be about 0.35m, which shows that the spatial resolution of the satellite image has an important influence on the DTM precision, and the higher the spatial resolution of the satellite image is, the higher the DTM construction precision and reliability are.

FIG. 9 also illustrates that in Chongming North beach, the time span also has some effect on the tidal beach DTM. DTM_LIDAnd annual DTMON-fusion (i.e., DTM)_LIDDTMnon-fusion) and statistical analysis of pixel-by-pixel differences are shown in fig. 9. In FIG. 9(b), 12% of the pixels are distributed in the range of-0.10 m-0.10m, 76% are distributed in the range of-0.30 m-0.30m, and RMSE is 0.27 m. Comparing the error statistics of fig. 6(c) and fig. 9(b), the absolute error of 24% of the pixels in DTMnon-fusion exceeds 0.3m, and the absolute error of 14% of the pixels in dtmffusion exceeds 0.3m, which indicates that the error of DTMnon-fusion includes the offset caused by the elevation change of the terrain itself within one year.

Claims

1. a single-day tidal flat digital terrain construction method of high time-space satellite data fusion, is characterized in that, this method comprises the following concrete steps:

Step 1, high spatiotemporal resolution image fusion, using the spatiotemporal adaptive reflection fusion model STARFM to generate hourly images with a spatial resolution of 16m;

Step 2, the tidal dynamic model TPXO8 obtains the tidal level, and the tidal dynamic model TPXO8 is used to obtain the tidal data;

Step 3, based on the extraction of the water edge of the BP neural network, using the BP neural network supervised classification method to accurately classify the muddy tidal wetlands and water bodies, and then extract the water edge;

Step 4, DTM construction, combined with tide level and waterline data to construct tidal flat terrain elevation.

2. a kind of single-day tidal flat digital terrain construction method of high time-space satellite data fusion according to claim 1, is characterized in that, described step 1 specifically comprises:

Step 1.1, first introduce the information of adjacent and similar pixels in the window, and use the _STARFM algorithm to first resample the GOCI image after geometric correction and georeferencing to the GF-1/WFV image with a resolution of 16m. , y _j ) on the type of ground cover and systematic errors are constant at prediction time t ₀ and image date t _k , then:

GF( _xi , y _j , t ₀ )=G( _xi , y _j , t ₀ )+GF( _xi , y _j , t _k )-G( _xi , y _j , t _k ) (1)

where GF and G represent GF-1/WFV and GOCI images, respectively, (x _i , y _j ) is the pixel position, t ₀ is the prediction time, and t _k is the time of the existing image pair;

Step 1.2, use the sliding window and weight function to predict the center pixel value

The STARFM algorithm assumes that for the same ground target, the error of GF-1/WFV and GOCI images is constant, and the error is systematic in a short time and only depends on the characteristics of a single pixel; therefore, using a certain time A pair of GF-1/WFV-GOCI images is used to estimate the error of each pixel, and then the GOCI images at other times are used to predict the high spatial resolution images of GF-1/WFV;

The STARFM algorithm first resamples the GOCI image after geometric correction and geo-referencing to the GF-1/WFV image with a resolution of 16m. Set the ground coverage type and system error on the pixel (x _i , y _j ) at the prediction time t ₀ and the image date t _k unchanged, then:

However, GOCI images are usually mixed pixels at GF-1/WFV 16m resolution, or within the forecast date, the tidal flat waterline is often changed due to tidal effects; by introducing additional information of adjacent pixels, the STARFM algorithm is used. The weighting function in calculates the surface reflectance of the pixel at the center of predicted time t ₀ :

where ω is the size of the search window, (x _ω/2 , y _ω/2 ) is the center pixel position of the search window; the weight W _ijk determines the contribution of each adjacent pixel to the estimated reflectivity of the center pixel; W _ijk is determined by the spectral difference, the time difference, the distance between the central pixel position (x _ω/2 , y _ω/2 ) and the candidate pixel position ( _xi , y _j );

The spectral differences between the GF-1/WFV and GOCI data at a given location are:

S _ijk =|GF( _xi , y _j , t _k )-G( _xi , y _j , t _k )| (3)

Equation (3) approximates the difference between the GF-1/WFV observations and the average value of adjacent pixels at low spatial resolution; a smaller value of S _ijk indicates that the pixel has a value closer to the average value of its surrounding pixels. Spectral features, i.e. the variation of high-resolution pixels is comparable to the average variation of surrounding pixels, at this time, a higher weight is assigned to the reflectance of the pixel in Eq. (2); the difference between the input GOCI data and the predicted GOCI data is The time difference is:

T _ijk =|G( _xi , y _j , t _k )-G( _xi , y _j , t ₀ )| (4)

T _ijk measures the change between prediction and acquisition time, a smaller T _ijk means less change between time t _k and t ₀ , pixels should be assigned higher weights; center pixel position at date t _k The distance between (x _ω/2 , y _ω/2 ) and the candidate pixel position ( _xi , y _j ) is:

d _ijk is the spatial distance between the central predicted pixel and the surrounding spectrally similar candidate pixels, and for closer pixels, assign a higher weight to its candidate pixels.

3. a kind of single-day tidal flat digital terrain construction method of high time-space satellite data fusion according to claim 1, is characterized in that, described step 2 specifically comprises:

In step 2.1, the tidal heights of the TPXO8 tidal model are given as complex amplitudes, so the subtidal division of a single component of frequency f at location x at time t is given by:

h(t,x)=pu(t,x)·Re[h(x)exp{i[f(tt ₀ )+V ₀ (t ₀ )+ph(t,x)]}] (6)

where V ₀ (t ₀ ) is the astronomical parameter at time t ₀ , and pu(t, x) and ph(t, x) are the node corrections; then the amplitude is |h|, and the phase is tan ^-1 (-Im(h )/Re(h)), IM and RE represent the imaginary number of the complex number, namely the Imagic Number part and the real number, namely the RealNumber part, i*i=-1;

Step 2.2, get the Matlab toolbox code from OSU-TPXO Tide Models (https://www, tpxo.net/home), and input the precise position coordinates and time information into the TPXO8 model to get the tide level data.

4. the single-day digital terrain construction method of a kind of high time-space satellite data fusion according to claim 1, is characterized in that, described step 3 specifically comprises:

Step 3.1, select two types of interest areas of intertidal wetlands and water bodies as training samples, so that the training samples are used as the input layer of supervised classification;

Step 3.2, use the BP neural network based on ENVI5.5 (trial version) to supervise and classify the images;

Step 3.3, output the supervised classification results and the degree of separation of the samples. Usually, if the difference between the two classification areas is obvious, the greater the degree of divisibility of the sample, and the better the degree of separation, the maximum degree of separation is 2.0, and the two classification areas The dividing line between them is the water edge.

5. the single-day digital terrain construction method of a kind of high time-space satellite data fusion according to claim 1, is characterized in that, described step 4 specifically comprises:

Step 4.1, construct Delaunay irregular triangulation, and use Delaunay irregular triangulation to divide;

Step 4.2, convert the TIN to grid, and interpolate the value into the grid DTM to obtain the digital terrain elevation.