CN102831310B - Method for building high-spatial resolution NDVI (normalized difference vegetation index) time series data - Google Patents
Method for building high-spatial resolution NDVI (normalized difference vegetation index) time series data Download PDFInfo
- Publication number
- CN102831310B CN102831310B CN201210295918.4A CN201210295918A CN102831310B CN 102831310 B CN102831310 B CN 102831310B CN 201210295918 A CN201210295918 A CN 201210295918A CN 102831310 B CN102831310 B CN 102831310B
- Authority
- CN
- China
- Prior art keywords
- ndvi
- modis
- pixel
- pixels
- time
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Expired - Fee Related
Links
- 238000000034 method Methods 0.000 title claims abstract description 36
- 238000005259 measurement Methods 0.000 claims description 9
- 238000004364 calculation method Methods 0.000 claims description 5
- 238000011049 filling Methods 0.000 claims description 2
- 238000003384 imaging method Methods 0.000 abstract description 2
- 230000008859 change Effects 0.000 description 10
- 230000007613 environmental effect Effects 0.000 description 8
- 230000002123 temporal effect Effects 0.000 description 7
- 230000000694 effects Effects 0.000 description 6
- 230000004927 fusion Effects 0.000 description 6
- 238000002310 reflectometry Methods 0.000 description 6
- 238000012544 monitoring process Methods 0.000 description 5
- 238000010586 diagram Methods 0.000 description 4
- 238000011160 research Methods 0.000 description 4
- 238000004458 analytical method Methods 0.000 description 3
- 238000002474 experimental method Methods 0.000 description 3
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 3
- 230000003044 adaptive effect Effects 0.000 description 2
- 238000000354 decomposition reaction Methods 0.000 description 2
- 238000001514 detection method Methods 0.000 description 2
- 230000003203 everyday effect Effects 0.000 description 2
- 230000007774 longterm Effects 0.000 description 2
- 239000000203 mixture Substances 0.000 description 2
- 238000007500 overflow downdraw method Methods 0.000 description 2
- 230000008569 process Effects 0.000 description 2
- 230000005855 radiation Effects 0.000 description 2
- 239000002689 soil Substances 0.000 description 2
- 239000002028 Biomass Substances 0.000 description 1
- 241000196324 Embryophyta Species 0.000 description 1
- 241000735234 Ligustrum Species 0.000 description 1
- 241001444165 Scolopax rusticola Species 0.000 description 1
- 241000209140 Triticum Species 0.000 description 1
- 235000021307 Triticum Nutrition 0.000 description 1
- 238000010521 absorption reaction Methods 0.000 description 1
- 230000004075 alteration Effects 0.000 description 1
- 235000013405 beer Nutrition 0.000 description 1
- 230000008901 benefit Effects 0.000 description 1
- 230000015572 biosynthetic process Effects 0.000 description 1
- 230000002354 daily effect Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 238000013507 mapping Methods 0.000 description 1
- 239000011159 matrix material Substances 0.000 description 1
- 238000002156 mixing Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000029553 photosynthesis Effects 0.000 description 1
- 238000010672 photosynthesis Methods 0.000 description 1
- 238000001556 precipitation Methods 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 238000012552 review Methods 0.000 description 1
- 239000011435 rock Substances 0.000 description 1
- 230000003595 spectral effect Effects 0.000 description 1
- 238000003786 synthesis reaction Methods 0.000 description 1
- 230000005068 transpiration Effects 0.000 description 1
- 238000012795 verification Methods 0.000 description 1
Landscapes
- Image Processing (AREA)
Abstract
The invention discloses a method for building high-spatial resolution NDVI (normalized difference vegetation index) time series data. The high-spatial resolution NDVI time series data is predicted and built according to low-spatial resolution MODIS (moderate-resolution imaging spectroradiometer) pixels in known MODIS NDVI data and high-spatial resolution TM (thematic mapper) pixels in TMNDVI (thematic mapper normalized difference vegetation index) data. TM data is combined with MODIS data, and accordingly high-spatial resolution NDVI time series data with quite fine precision can be obtained effectively.
Description
Technical Field
The present invention relates to a method for constructing high spatial resolution NDVI time series data, which can generate high spatial resolution NDVI time series data through known MODIS NDVI data and TM NDVI data.
Background
The Normalized difference vegetation index (abbreviated as NDVI) is a widely used vegetation index, which can be obtained from the earth surface reflectivity of the red and near infrared bands obtained by satellite remote sensing (see appendix 17), for example, TM and MODIS data of the united states terrestrial satellite include NDVI data obtained by satellite remote sensing, which can be downloaded or purchased from related websites.
NDVI can be used to detect vegetation growth status, vegetation coverage, etc., which can reflect the background effects of plant canopy, such as soil, wet ground, snow, etc. parameters related to vegetation coverage. NDVI is more than or equal to-1 and less than or equal to +1, and a negative value indicates that the earth surface or an aerial region corresponding to the earth surface is cloud, water, snow and the like, and has high reflection to visible light; 0 represents rock, bare earth, or the like; positive values indicate vegetation coverage and increase with increasing coverage.
The NDVI data obtained by satellite detection is placed on a map to form a digitized NDVI value picture, each pixel point in the picture corresponds to a square area of the ground surface, the square area has a specific NDVI value, all the pixel points are arranged together by the NDVI values of the pixel points to form the digital picture expressed by the NDVI values, and the NDVI value of each pixel of the picture can be used for expressing the ground surface vegetation condition of the position on the map.
The TM NDVI data refers to surface data obtained by a Thermal Mapper (TM) sensor mounted on a landdetecting satellite system Landsat developed by the national aerospace agency, and each surface pixel in the TM NDVI data has a corresponding NDVI value. Similarly, the MODIS NDVI data refers to surface data obtained by a medium-resolution imaging spectrometer (abbreviated as MODIS) mounted on the Aqua and Terra of the american earth observation system series, where each surface pixel in the MODIS NDVI data has a corresponding NDVI value.
The spatial resolution of the NDVI value in the TM NDVI data is 30m × 30m, and the temporal resolution is 16 days, while the maximum spatial resolution of the NDVI value in the MODIS NDVI data is 250m × 250m, and the temporal resolution is 1 day (the MODIS NDVI data may actually be provided by two satellites, the spatial resolution and the temporal resolution of the MODIS NDVI data obtained by the two satellites are the same, except that the two satellites are respectively in the morning and afternoon, so that for the same area, the MODIS NDVI data may be obtained twice a day, and for the sake of simplicity, the MODIS NDVI data in the same area are selected in the same area in the same day of the morning. That is, for the same location on the earth, such as somewhere on the loess plateau, the Landsat satellite passes every 16 days, and each pixel of the surface image obtained by the TM sensor thereon represents a scale range of 30m × 30m, so that the spatial resolution of NDVI data obtained by the TM sensor on a spatial scale is 30m × 30 m. Similarly, the Aqua or Terra satellites will pass through the site once a day, and each pixel of the image obtained by the MODIS sensor above it represents a scale range of 250m by 250m, so that the spatial resolution of the NDVI data obtained by the MODIS sensor on the spatial scale is 250m by 250 m.
As can be seen from the above description of the TM NDVI data and the MODIS NDVI data, the TM NDVI data has a higher spatial resolution and a lower temporal resolution, while the MODIS NDVI data has a lower spatial resolution and a higher temporal resolution. The spatial resolution of TM NDVI data is higher, the detection precision of the TM NDVI data to the ground vegetation situation is much better than that of MODIS NDVI data, but the defects are that the provided NDVI values are not available every day, the NDVI values can be acquired only once at intervals of 16 days, the continuous monitoring of a specific area has difficulty, and the vegetation change situation in an accurate time range cannot be acquired, for example, the monitoring of forest fires in a certain mountain area or the monitoring of geological disasters has a blind spot in time. On the other hand, although the MODIS NDVI data can be acquired every day, the spatial resolution is not sufficient, and the data of the changing situation cannot be provided for the forest fire or geological disaster in a small range.
Therefore, it is desirable to provide a method for constructing high-spatial-resolution NDVI time-series data by predicting known MODIS pixels with low spatial resolution in the MODIS NDVI data and TM pixels with high spatial resolution in the TM NDVI data
Disclosure of Invention
In order to solve the technical problem, the invention provides a method for constructing high-spatial-resolution NDVI time series data, which is used for predicting and constructing the high-spatial-resolution NDVI time series data through known MODIS pixels with low spatial resolution in the MODIS NDVI data and TM pixels with high spatial resolution in the TM NDVI data.
The method comprises the following steps:
(A) will t1Time t and2the temporally low spatial resolution MODIS pixels are filled at the TM pixel scale and thenDividing the TM pixels into a feature type c with a total type of l according to an NDVI threshold value, and setting that the TM pixels of the same feature type c have the same NDVI increase rate; the pixels of each "picture" in the TM NDVI data are referred to as the TM pixels; the pixels of each "picture" in the MODIS NDVI data are referred to as the MODIS pixels;
(B) measurement obtained from t1To t2The NDVI growth rate of the MODIS pixel at a time;
(C) adding the NDVI growth rates of the TM pixels of different surface feature types in the surface feature type with the total class l to the area ratio of the TM pixels of the surface feature type in the MODIS pixel, so as to obtain the average value of the NDVI growth rates of all the TM pixels filled in the MODIS pixel, wherein the average value is equal to the NDVI growth rate of the MODIS pixel in the step (B);
(D) setting the NDVI increasing rate of the same ground object type in at least (l-1) MODIS pixels adjacent to the MODIS pixel, and repeating the steps (A) to (C) for the (l-1) MODIS pixels so as to calculate and obtain the NDVI increasing rate of the TM pixel of each ground object type;
(E) measurement to obtain t1Setting the NDVI value of the TM pixel at the moment to be linearly changed along with time, and calculating to obtain t according to the NDVI increasing rate of the TM pixel of each ground feature type obtained in the step (D)2NDVI value of TM pixel at time;
(F) t obtained by the calculation2NDVI value of TM pixel at time instant with said known t1The NDVI values of the TM pixels at a time are temporally aligned, and the high spatial resolution NDVI time series data is obtained.
Preferably, in step E, the set NDVI value varies linearly with time according to formula (1):
NDVI2＝NDVI1+k×(t2—t1) (1)
wherein NDVIiIs ti(i ═ 1,2) the NDVI value at time; k represents a number from t1To t2The NDVI growth rate corresponding to the time instant.
Preferably, in step C, the average value of the NDVI increase rates of all TM pixels is equal to the NDVI increase rate of the MODIS pixel, which is obtained from formula (2):
kMODIS(x,y,t1→t2)＝∑l c＝1fc(x,y,t1)×kc TM(x,y,t1→t2) (2)
wherein k isMODIS(x,y,t1→t2) Is from t1To t2The growth rate of the MODIS pixel (x, y) at time; k is a radical ofc TM(x,y,t1→t2) Is from t1To t2The growth rate of the corresponding c-th type TM pixel in the time MODIS pixel (x, y); f. ofc(x,y,t1) Is t1The area ratio of the class c feature in the time pixel (x, y); l is the total ground object type in pixel (x, y).
Preferably, in step D, the NDVI growth rate of the TM pixel of each surface feature class is obtained by formula (3):
wherein Δ t represents the value from t1To t2The time of day.
Preferably, the ground object class i is 5 or 6.
The invention provides a method for constructing high-spatial-resolution NDVI time sequence data, which predicts a high-spatial-resolution NDVI numerical value according to the relevance of MODIS NDVI data and TM NDVI data, so that the high-spatial-resolution NDVI time sequence data can be constructed. The method integrates TM data and MODIS data, can efficiently obtain NDVI time sequence data with higher resolution and has quite good precision.
Drawings
The drawings are only for purposes of illustrating and explaining the present invention and are not to be construed as limiting the scope of the present invention. Wherein,
FIG. 1 schematically shows the correspondence between TM NDVI data and MODIS NDVI data;
FIG. 2 is a schematic diagram of the MODIS NDVI data prediction based on the MODIS NDVI data in FIG. 1;
fig. 3 is a diagram illustrating a correspondence relationship between MODIS pixels and TM pixels.
Detailed Description
In order to more clearly understand the technical features, objects, and effects of the present invention, embodiments of the present invention will now be described with reference to the accompanying drawings. Wherein like parts are given like reference numerals.
Fig. 1 schematically shows the corresponding situation of TM NDVI data and MODIS NDVI data, and in fig. 1, the TM NDVI data and the MODIS NDVI data are visually represented in the form of "pictures" to represent the difference between the temporal resolution and the spatial resolution.
It will be understood by those skilled in the art that, as described above, the TM NDVI data and the modistvi data in fig. 1 are imaginable in a "picture" format, wherein the TM NDVI data is formed by a group of digital "pictures" arranged in a time sequence, wherein adjacent "pictures" are taken at a time interval of 16 days, and each pixel in each "picture" corresponds to a 30m by 30m square region of the surface, which has its specific NDVI value; similarly, MODIS NDVI data is also made up of a set of digital "pictures" arranged in a chronological order, where adjacent "pictures" are taken at 1 day intervals, and each pixel in each "picture" corresponds to a 250m by 250m square region of the earth's surface with its particular NDVI value.
As shown, TM NDVI data is shown above the time axis, that is, a "picture" with a spatial resolution of 30m x 30m is obtained every 16 days with a certain time point as an origin, and thus, the pixels of each "picture" in the TM NDVI data can be referred to as TM pixels, each TM pixel having an NDVI value corresponding thereto. While below the time axis, MODIS NDVI data is shown, with the same time point as the origin (in practice, MODIS NDVI data is less than 30 minutes different from the origin, but the delay is negligible for 24 hours per day), a "picture" with a spatial resolution of 250m x 250m is obtained every 1 day interval, and therefore, the pixels of each "picture" in MODIS NDVI data can be referred to as MODIS pixels, each having an NDVI value corresponding thereto.
The invention provides a method for constructing high-spatial-resolution NDVI time sequence data, which is characterized in that high-spatial-resolution NDVI values of TM NDVI data at the same moment are predicted according to the relevance of MODIS NDVI data and TM NDVI data and low-spatial-resolution NDVI values of MODIS NDVI data, and finally the predicted high-spatial-resolution NDVI values of the TM NDVI data are arranged according to time, so that the high-spatial-resolution NDVI time sequence data can be constructed.
The concrete description is as follows:
example 1
Referring to fig. 2, a schematic diagram of the prediction of the TM NDVI data according to the MODIS NDVI data based on fig. 1 is shown, wherein the TM NDVI data obtained by prediction is shown by a dotted line.
Since the satellites from which the TM NDVI data and the MODIS NDVI data are obtained have similar orbital parameters and a satellite transit time interval of less than 30 minutes, it can be theoretically assumed that the TM NDVI data and the MODIS NDVI data are linearly variable in a short time in time series.
That is, in fig. 2, the NDVI values of TM pixels adjacent to "pictures" in the TM NDVI data above the time axis are assumed to change linearly with time, and since there is similarity between the MODIS NDVI data below the time axis and the TM NDVI data from the satellite orbit to the transit time, it can be similarly inferred that the ndi values of MODIS pixels adjacent to "pictures" in the MODIS NDVI data below the time axis also change linearly.
Thus, for two different times t1、t2The linear change relationship between the NDVI values can be expressed as formula (1):
NDVI2＝NDVI1+k×(t2—t1) (1)
wherein NDVIiIs ti(i ═ 1,2) the NDVI value at time; k represents a number from t1To t2The NDVI growth rate corresponding to the time instant.
That is, the above formula (1) is common to both MODIS NDVI data and TM NDVI data, and should have an associated growth rate k for both. Therefore, if only the NDVI growth rate k of the MODIS pixel is obtained, the NDVI growth rate of the TM pixel is obtained, and the NDVI growth rate is substituted into the formula (1) to predict t2NDVI value of TM pixel at time instant.
That is, we can convert the NDVI increase rate k of 1-day-spaced MODIS pixels into the NDVI increase rate of TM pixels, and can also calculate the NDVI values of 1-day-spaced TM pixels, so that the NDVI values of high spatial resolution can also be 1-day-spaced, and as shown by the blocks shown by the dotted lines in fig. 2, the calculated NDVI values of high spatial resolution (NDVI values of TM pixels) are represented, and the high spatial resolution NDVI time series data are obtained by time arrangement.
In particular to implementation, because the MODIS NDVI data and the TM NDVI data have different spatial resolutions, the NDVI growth rate of the MODIS pixels shows the NDVI growth condition of a 250m by 250m area on the ground, and the NDVI growth rate of the TM pixels shows the NDVI growth condition of a 30m by 30m area on the ground. Therefore, it is conceivable that only 1 MODIS pixel is contained in the same land area, for example, the land area 250m × 250m, but approximately 8 × 64 TM pixels are contained, so that the NDVI increase rate of 1 MODIS pixel in the area should be substantially equal to the average value of the NDVI increase rates of 64 TM pixels in the area.
Thus, the present invention provides a solution to first fill the low spatial resolution MODIS pixels at the TM pixel scale, e.g., t below the time axis in FIG. 21Time t and2the MODIS pixels at a time are filled according to a TM pixel scale, then the TM pixels are divided into different categories, the TM pixels in the same category have the same NDVI growth rate, and the increase rate is measured from t1To t2The NDVI increase rates of TM pixels at the same position in the MODIS pixels adjacent to each other at the time are added after the measured NDVI increase rates of the TM pixels of different classes are multiplied by the area ratio occupied by the TM pixel of the class in the MODIS pixel where the TM pixel is located, so as to obtain the average value of the NDVI increase rates of all TM pixels filled in the MODIS pixel, and the average value should be equal to the NDVI increase rate of the MODIS pixel.
That is, the relationship between the NDVI value of the MODIS pixel and the NDVI increase rate of the NDVI value of the TM pixel can be expressed by equation (2):
kMODIS(x,y,t1→t2)＝∑l c＝1fc(x,y,t1)×kc TM(x,y,t1→t2) (2)
wherein k isMODIS(x,y,t1→t2) Is from t1To t2The growth rate of the MODIS pixel (x, y) at time instant, which can be obtained by measurement; k is a radical ofc TM(x,y,t1→t2) Is from t1To t2The growth rate of the corresponding c-th type TM pixel in the time MODIS pixel (x, y), which is the unknown parameter to be obtained; f. ofc(x,y,t1) Is t1Area of class c feature in time pixel (x, y)This can be obtained by measurement; l is the total ground object type in the pixel (x, y), which is a value set according to the surface condition.
That is, we assume that pixels with similar NDVI have the same growth rate, so that TM pixels can be divided into l terrain types, for example, TM pixels can be divided into c different terrain types according to the difference of thresholds of different terrain expressing NDVI values, considering that the NDVI of the water body is lower than 0, the value of the NDVI of the bare terrain ranges from 0 to 0.2, the NDVI of the low-coverage vegetation ranges from 0.2 to 0.5, and the NDVI of the high-coverage vegetation is substantially above 0.5, so the threshold classification can be used as follows: water at 1.0-0, bare soil at 0-0.2, low vegetation coverage at 0.2-0.5, and dense vegetation at 0.5-1.0, to avoid bias from a single ground type.
The formula (2) makes it possible to obtain the growth rate of the TM pixel using the MODIS pixel. However, there is still a problem that the total number of surface feature types l is often greater than 1, which means that the value of the unknown parameter cannot be obtained unless there are at least another (l-1) equations.
To solve this problem, MODIS neighboring pixels are introduced to provide more information, assuming that the growth rates of the same type of feature are consistent among the neighboring pixels. The linear mixture model of N neighboring pixels has formula (3) in the form of a generic matrix:
wherein Δ t represents the value from t1To t2Time of day, k in equation (3)MODIS(x, y,. DELTA.t) and fc(x, y) can be obtained according to the existing MODIS time series data and TM single-time data. Thus, kc TM(x, y,. DELTA.t) can be estimated by solving a system of linear equations. Finally, according to the formula (1), the increase rate of the TM pixel is used for calculating and obtaining t2The NDVI value at that time is predicted to be the NDVI value shown by the dashed line in fig. 2.
The above description is given by way of example in a manner of intuitively adjoining pixels of an NDVI digital "picture" according to the characteristics of TM NDVI data and MODIS NDVI data, and therefore, for the entire NDVI digital "picture", as long as the above steps are traversed for each MODIS pixel and TM pixel, each pixel for obtaining TM NDVI data can be calculated at t2NDVI value at time.
Similarly, as shown in FIG. 2, for t2TM pixels after a time instant, e.g. t3The TM pixel at time can be adjacent t by the same method2To t3The NDVI increase rate of the pixels of the time MODIS is based on the predicted t2NDVI value of TM pixel at time is further calculated to obtain t3NDVI value of TM pixel at time, and so on. Finally, t obtained by prediction is sequentially obtained1、t2、t3… …, etc., are arranged in time to obtain a series of high spatial resolution NDVI time series data.
Based on the above analysis, the method of the present invention can be summarized as including the steps of:
(A) will t1Time t and2filling MODIS pixels with low spatial resolution at a moment according to TM pixel scale, dividing the TM pixels into a feature type c with a total type of l according to an NDVI threshold value, and setting the TM pixels of the same feature type c to have the same NDVI increase rate;
(B) measurement obtained from t1To t2The NDVI growth rate of the MODIS pixel at a time;
(C) adding the NDVI increase rates of the TM pixels of the different feature types to the area ratio of the TM pixel of the feature type in the MODIS pixel where the TM pixel is located, so as to obtain an average value of the NDVI increase rates of all the TM pixels filled in the MODIS pixel, where the average value is equal to the NDVI increase rate of the MODIS pixel in the step (B);
(D) setting the NDVI increasing rate of the same ground object type in at least (l-1) MODIS pixels adjacent to the MODIS pixel, and repeating the steps (A) to (C) for the (l-1) MODIS pixels so as to calculate and obtain the NDVI increasing rate of the TM pixel of each ground object type;
(E) measurement to obtain t1Setting the NDVI value of the TM pixel at the moment to be linearly changed along with time, and calculating to obtain t according to the NDVI increasing rate of the TM pixel of each ground feature type obtained in the step (D)2NDVI value of TM pixel at time;
(F) t obtained by the calculation2NDVI value of TM pixel at time instant with said known t1The NDVI values of the TM pixels at a time are temporally aligned, and the high spatial resolution NDVI time series data is obtained.
In step E, the set NDVI value varies linearly with time as obtained from equation (1).
Wherein in step C, the average value of the NDVI increase rates of all TM pixels is equal to the NDVI increase rate of the MODIS pixel, which is obtained from formula (2).
Wherein in step D, the NDVI growth rate of the TM pixel of each surface feature class is obtained by formula (3).
Example 2
The method of constructing high spatial resolution NDVI time series data according to the present invention is further described below in one specific example. As shown in fig. 3, a diagram of correspondence between MODIS pixels and TM pixels is shown.
Referring to fig. 3, the left side of the figure shows a 3 × 3 MODIS pixel, and each small square indicates a MODIS pixel, which is 9 MODIS pixels in total. The right side of the figure shows an enlarged view of the left side labeled MODIS pixels (target MODIS pixels) filled on the TM pixel scale, i.e., each small square on the right side represents a corresponding one of the TM pixels. It should be noted that, since the spatial resolution of the MODIS pixels is 250m × 250m and the spatial resolution of the TM pixels is 30m × 30m, for convenience of later calculation, the MODIS pixels are resampled to 240m × 240m so that each MODIS pixel is equal to 8 × 8 TM pixels, and each MODIS pixel is divided into 64 TM pixels by 8 in the drawing for simplicity of description, that is, the target MODIS pixel on the left side in the drawing is filled with 64 TM pixels according to the TM pixel scale.
First, according to a known time t1The NDVI values of the TM pixels of (1) are classified by using threshold values 0 and 0.2, where the NDVI is less than 0 and is a first type of ground feature, 0-0.2 is a second type of ground feature, and more than 0.2 is a third type of ground feature, and the obtained classification results are shown in the figure, where the target MODIS pixel includes 24 TM pixels of the 1 st type of ground feature, 20 TM pixels of the 2 nd type of ground feature, and 20 TM pixels of the 3 rd type of ground feature.
According to the classification, the area proportion of each type of ground object in the MODIS pixels can be calculated. According to FIG. 3, the area ratios of the three types of ground objects in the target MODIS pixel are calculated to be f1＝0.375，f2＝0.3125，f3＝0.3125。
According to formula (2):
kMODIS(x,y,t1→t2)＝∑l c＝1fc(x,y,t1)×kc TM(x,y,t1→t2) (2)
for the present example, the simplified representation of the pixel from t1To t2The rate of change at a time can be obtained as follows:
kMODIS＝f1*k1 TM+f2*k2 TM+f3*k3 TM (4)
in equation (4), k is due to the left sideMODISCan be measured by measuring t1Time t2The NDVI growth rate of the MODIS pixel at a time is obtained, so it is necessary to obtain three unknowns k on the right side of the equation1 TM、k2 TM、k3 TMHowever, one MODIS pixel can only obtain the above-mentioned one equation (4) related to the three unknowns, so in order to obtain the three unknowns, more equations need to be introduced.
According to the assumption that the NDVI growth rates of the same kind of feature in the neighboring pixels are consistent, as shown in fig. 3,2 pixels closest to the target MODIS pixel (e.g., a pixel a above and a pixel c below the target MODIS pixel b on the left side in the example drawing) are selected. Then, the TM pixel NDVI change rate from the observation time (t1) to the prediction time (t2) is calculated pixel by pixel.
According to the NDVI change rate of the TM pixels obtained in the above steps and the ratio of each category corresponding to each TM pixel, 3 equations can be listed as follows:
ka MODIS＝fa,1*k1 TM+fa,2*k1 TM+fa,3*k3 TM (4a)
kb MODIS＝fb,1*k1 TM+fb,2*k2 TM+fb,3*k3 TM (4b)
kc MODIS＝fc,1*k1 TM+fc,2*k2 TM+fc,3*k3 TM (4c)
the pixel a above the target MODIS pixel corresponds to equation (4a), the target MODIS pixel b itself corresponds to equation (4b), and the pixel c below the target MODIS pixel corresponds to equation (4 c).
The rate of change of the MODIS pixels to the left of the above three equations is known, t1The area ratio of each type of ground object at the moment is also known, so that three unknowns k can be obtained by solving the three equations1 TM、k2 TM、k3 TM。
In the actual operation process, in order to improve the calculation accuracy, an equation can be selected for each of the 8 MODIS pixels a, c, d1, d2, d3, e1, e2, and e3 closest to the target pixel, so that 9 equations can be obtained, and as long as the feature type l is less than or equal to 9, the NDVI change rate corresponding to each feature type can be obtained through a common least square method.
And substituting the NDVI change rate corresponding to each ground feature type obtained by solving into a formula (1), so as to obtain the NDVI value of each ground feature type at the time t 2. The values are temporally aligned to obtain high spatial resolution NDVI time series data.
The application provides a method for generating high-spatial-resolution NDVI time sequence data based on MODIS NDVI time sequence data and TM NDVI data, the method combines an NDVI linear growth model and an expanded NDVI linear mixed model together, and converts an original prediction problem into a linear system solving problem. Since in a linear system it is estimated by several parameters (the growth rate of various types of TM pixels), the method introduces available information from neighboring pixels based on the assumption that the same class growth rate of neighboring pixels is consistent.
This method has proven to be efficient and accurate according to several experiments we performed. Firstly, the experimental result shows that the influence of the category number on the result comprises two aspects of precision and texture information, the two aspects are considered in a comprehensive mode, and 5-6 categories are recommended to be used for prediction based on the existing research area data.
Secondly, there may be some difference in the performance of the method for homogeneous and heterogeneous regions. We performed several experiments on two different regions, one of which was more homogeneous and the other was more heterogeneous. Experimental results for the homogeneous region show that the accuracy of the predictions of the present invention is at least comparable to, and even many times higher than, the accuracy of the ESTARFM method (see appendix 16). At the same time, experimental results in the heterogeneous region show that the new method can give better results than ESTARFM.
In addition, in this method, TM pixels need to be classified based on NDVI values, and there is a potential problem of a phase effect in the classification. Classification based on data from different time phases results in different predictions, and therefore we performed several experimental verifications of the problem. The results of these experiments indicate that the temporal effects have a significant impact in long-term predictions, and that appropriate images must be selected to minimize errors, so that it is recommended that TM NDVI data be used at the time of best vegetation growth when long-term predictions are made. However, for a short period of time, the effect of the phase effect is of secondary importance, since the pixels with similar NDVI values vary similarly in a short period of time, and there is no specific phase requirement for the TMNDVI data.
The method proposed by the present application also has a very important advantage, namely efficiency. The ESTARFM method (see appendix 16) has good accuracy in fusion, but its time cost is extremely high, taking 2-3 hours to process a phase of TM data, and at least 40 hours if annual data is desired, making its application impractical on a global or regional scale. The method provided by the application only needs 1.5-2 hours for processing the annual data of the TM data, and the high-efficiency property enables the method to be suitable for being popularized and applied to ecological models on a global or regional scale.
The above description is only an exemplary embodiment of the present invention, and is not intended to limit the scope of the present invention. Any equivalent alterations, modifications and combinations can be made by those skilled in the art without departing from the spirit and principles of the invention.
Appendix of the prior art:
[1] the spectral reflectance of wheat was used to estimate the absorption and radiation as well as the leaf area index, Journal of agriculture (1984), 76, 300-.
[2] The combination of medium and low spatial resolution satellite data improves the sub-pixel NDVI time series estimation, environmental Remote Sensing, 112(2008), 118-.
[3] A simple method for reconstructing a high quality NDVI time series data set based on a Savitzky-Golay filter, environmental Remote Sensing, 91(2004) 332-.
[4] Inversion of reflectivity time profiles in heterogeneous landscape areas using simulated and real NOAA-AVHRR data, Remote Sensing International Journal (2000),21,753-
[5] Estimate sub-pixel surface variables for the input parameters of the climate hydrological model, national center for space research (france) INRS-Eau research report (2000), vol.564.64pp (Fortin, j.p., Bernier, m., El Battay, a., & Gauthier, y.rapport de receherche INRS-Eau, vol.564.64pp).
[6] Fusing Landsat and MODIS surface reflectivity: the daily Landsat surface reflectance (Gao, F., Masek, J., Schwaller, M., & Hall, F. (2006). IEEE Transactions on Geosilence and Remote Sensing,44, 2207-.
[7] A new data fusion method for high-space-time resolution forest mapping based on Landsat and MODIS NDVI data is disclosed, wherein environment Remote Sensing (2009)113, 1613-.
[8] The space-time adaptive reflectivity fusion model is used for fusing MODIS NDVI data to generate dense time series synthesis Landsat data, and Environment Remote Sensing (2009) (Hilker, T., Wulder, M.A., Coops, N.C., Sritz, N.s., White, J.C., Gao, F., et al (2009); Remote Sensing of environmental, 113, 1988-doping 1999).
[9] The method comprises the steps of real-time vegetation dynamic downscaling research of fusion of multi-temporal MODIS and Landsat NDVI data in a complex terrain area, and Environment Remote Sensing (2011) (Hwang, T, Song, C, Bolstad, V.P, Band, E.L., (2011) Remote Sensing of environmental 115(2011) 2499-2512).
[10] The Argentina Panpa region was based on a linear mixture theory NOAA-AVHRR NDVI decomposition and sub-pixel classification study (Kerdiles, H., & Grondona, M.O. (1995). International journal removal Sensing,16, 1303-.
[11] High and low spatial resolution NDVI data fusion methods for mediterranean environmental vegetation monitoring (Maselli, f., Gilabert, m.a., & Conese, C. (1998). Remote sensing of Environment,63, 208-.
[12] A variety of surface coverage type reflectivity data mixed pixel decomposition studies were performed using low spatial resolution satellite data (Oleson, k., Sarlin, s., Garrison, j., Smith, s., privet, j., & Emery, W. (1995). Remote Sensing of Environment,54, 98-112).
[13] Canopy reflectance, photosynthesis and transpiration (Sellers, P.J (1985). International journal of Remote Sensing,6, 1335-.
[14] Linear mixing theory and estimation of the proportion of the types of surface coverage (Settle, J.J., & Drake, N.A. (1993). International Journal of Remote Sensing,14, 1159-.
[15] Linear fusion of red and infrared bands for vegetation monitoring (Tucker, C.J. (1979) Remote Sensing of environmental, 8, 127- & 150).
[16] An enhanced spatio-temporal adaptive reflectivity fusion model of complex heterogeneity region (Zhu, X.L., J.Chen, F.Gao, X.H.Chen & J.G.Masek (2010); Remote Sensing of environmental, 114, 2610-2623).
[17] Review of radiation and biophysical performance of MODIS vegetation index (Huete, A., K.Didan, T.Miura, E.P.Rodriguez, X.Gao & L.G.Ferreira (2002); Remote Sensing of environmental, 83, 195-.
[18] Alpine meadow climates were monitored using NDVI time series data from AVHRR, vegettion and MODIS-compared to actual measurements (Fontana, f., c.rixen, t.jonas, g.aberegg & s.wunderle (2008). Sensors,8,2833- > 2853).
[19] Evaluation of spatiotemporal changes of vegetation phenology using NDVI time series data and fourier analysis: results of the Nanobian drought grassland environment (Wagenseil, H. & C. Saimi (2006). International journal of Remote Sensing,27, 3455-.
[20] Influence of pre-season temperature and precipitation on the eastern grassland spring phenology in Qinghai-Tibet plateau (Shen, M.G., Y.H.Tang, J.Chen, X.L.Zhu & Y.H.Zheng (2011). Agricutural and Forest technology, 151, 1711-.
[21] 1981-.
[22] 1981-.
[23]1982-1999 high latitude area of northern hemisphere normalized difference vegetation index and growth quartering trend (Tucker, C.J., D.A.Slayback, J.E.Pinzon, S.O.Los, R.B.Myneni & M.G.Taylor (2001). International Journal of biomethology, 45, 184-190).
[24]1981-1999 northern hemisphere vegetation activity changes inferred using satellite vegetation index data (Zhou, L.M., C.J.Tucker, R.K.Kaufmann, D.Slayback, N.V.Shabanov & R.B.Myneni (2001). Journal of geographic Research-actuators, 106, 20069-20083).
[25] Egyptian agricultural land status: polymorphic NDVI features obtained using Landsat TM (Lenney, m.p., c.e.woodcock, j.b.collins & h.hamdi (1996). RemoteSensing of Environment,56, 8-20).
[26] Landsat and MODIS NDVI data fusion products (Walker, j.j., k.m.de beers, r.h.wynne & f.gao (2012). RemoteSensing of Environment,117,381-393) were evaluated for forest climate analysis in arid areas.
[27] High spatial-temporal resolution NDVI data were generated and used for crop biomass estimation (Meng, j., x.du & b.wu (2011).
Claims (5)
1. A method for constructing high spatial resolution NDVI time sequence data is provided, which predicts and constructs high spatial resolution NDVI time sequence data through known MODIS pixels with low spatial resolution in the MODIS NDVI data and TM pixels with high spatial resolution in the TM NDVI data, wherein the pixel of each 'picture' in the TM NDVI data is called the TM pixel; the pixels of each "picture" in the modistvi data are referred to as the MODIS pixels;
characterized in that the method comprises the following steps:
(A) will t1Time t and2filling MODIS pixels with low spatial resolution at a moment according to TM pixel scale, dividing the TM pixels into a feature type c with a total type of l according to an NDVI threshold value, and setting the TM pixels of the same feature type c to have the same NDVI increase rate;
(B) measurement obtained from t1To t2The NDVI growth rate of the MODIS pixel at a time;
(C) adding the NDVI growth rates of the TM pixels of different surface feature types in the surface feature type with the total class l to the area ratio of the TM pixels of the surface feature type in the MODIS pixel, so as to obtain the average value of the NDVI growth rates of all the TM pixels filled in the MODIS pixel, wherein the average value is equal to the NDVI growth rate of the MODIS pixel in the step (B);
(D) setting the NDVI increasing rate of the same ground object type in at least (l-1) MODIS pixels adjacent to the MODIS pixel, and repeating the steps (A) to (C) for the (l-1) MODIS pixels so as to calculate and obtain the NDVI increasing rate of the TM pixel of each ground object type;
(E) measurement to obtain t1Setting the NDVI value of the TM pixel at the moment to be linearly changed along with time, and calculating to obtain t according to the NDVI increasing rate of the TM pixel of each ground feature type obtained in the step (D)2NDVI value of TM pixel at time;
(F) t obtained by the calculation2NDVI value of TM pixel at time instant with said known t1The NDVI values of the TM pixels at a time are temporally aligned, and the high spatial resolution NDVI time series data is obtained.
2. The method of claim 1, wherein in step E, the set NDVI value varies linearly with time according to equation (1):
NDVI2＝NDVI1+k×(t2—t1) (1)
wherein NDVIiIs ti(i ═ 1,2) the NDVI value at time; k represents a number from t1To t2The NDVI growth rate corresponding to the time instant.
3. The method according to claim 1 or 2, wherein in step C, the average value of the NDVI increase rates of all TM pixels is equal to the NDVI increase rate of the MODIS pixel, which is obtained from formula (2):
wherein k isMODIS(x,y,t1→t2) Is from t1To t2The growth rate of the MODIS pixel (x, y) at time; k is a radical ofc TM(x,y,t1→t2) Is from t1To t2Time MODIS pixel(x, y) growth rate of the corresponding class c TM pixel; f. ofc(x,y,t1) Is t1The area ratio of the class c feature in the time pixel (x, y); l is the total ground object type in pixel (x, y).
4. The method according to claim 1 or 2, wherein the terrain category/is 5 or 6.
5. The method according to claim 3, wherein the terrain category/is 5 or 6.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201210295918.4A CN102831310B (en) | 2012-08-17 | 2012-08-17 | Method for building high-spatial resolution NDVI (normalized difference vegetation index) time series data |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201210295918.4A CN102831310B (en) | 2012-08-17 | 2012-08-17 | Method for building high-spatial resolution NDVI (normalized difference vegetation index) time series data |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN102831310A CN102831310A (en) | 2012-12-19 |
| CN102831310B true CN102831310B (en) | 2015-04-01 |
Family
ID=47334443
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201210295918.4A Expired - Fee Related CN102831310B (en) | 2012-08-17 | 2012-08-17 | Method for building high-spatial resolution NDVI (normalized difference vegetation index) time series data |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN102831310B (en) |
Families Citing this family (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US9401030B2 (en) | 2014-04-25 | 2016-07-26 | Tazco Soil Service Co. | Image processing system for soil characterization |
| CN105279515A (en) * | 2015-09-14 | 2016-01-27 | 哈尔滨工业大学 | Multi-temporal remote sensing image feature extraction method based on difference temporal index features |
| CN107103584B (en) * | 2017-04-11 | 2019-06-14 | 北京师范大学 | A method of the production high-spatial and temporal resolution NDVI based on space-time weighting |
| CN107782700B (en) * | 2017-09-07 | 2019-07-09 | 北京师范大学 | A kind of AVHRR Reflectivity for Growing Season method for reconstructing, system and device |
| CN109685108B (en) * | 2018-11-23 | 2022-05-03 | 电子科技大学 | Method for generating high-space-time resolution NDVI long-time sequence |
| CN110147832B (en) * | 2019-05-08 | 2021-03-05 | 北京师范大学 | Method model for obtaining high space-time sequence EVI (error vector index) based on Bayesian theory |
| CN110147617A (en) * | 2019-05-22 | 2019-08-20 | 西南大学 | One kind carrying out extended method to GIMMS data based on MODIS data |
| CN110909821B (en) * | 2019-12-03 | 2020-07-28 | 中国农业科学院农业资源与农业区划研究所 | Method for carrying out high-space-time resolution vegetation index data fusion based on crop reference curve |
| CN112052720B (en) * | 2020-07-16 | 2021-04-20 | 贵州师范学院 | High-space-time normalization vegetation index NDVI fusion model based on histogram clustering |
Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102609726A (en) * | 2012-02-24 | 2012-07-25 | 中国科学院东北地理与农业生态研究所 | Method for classifying remote sensing images blended with high-space high-temporal-resolution data by object oriented technology |
Family Cites Families (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20120101784A1 (en) * | 2010-10-25 | 2012-04-26 | Trimble Navigation Limited | Wide-area agricultural monitoring and prediction |
-
2012
- 2012-08-17 CN CN201210295918.4A patent/CN102831310B/en not_active Expired - Fee Related
Patent Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102609726A (en) * | 2012-02-24 | 2012-07-25 | 中国科学院东北地理与农业生态研究所 | Method for classifying remote sensing images blended with high-space high-temporal-resolution data by object oriented technology |
Non-Patent Citations (3)
| Title |
|---|
| Estimation of planted forest leaf area index from TM imagery using the algorithm based on geometric-optical model;Hanyue Chen等;《2012 Second International Workshop on Earth Observation and Remote Sensing Applications》;20120611;第230-234页 * |
| 万华伟等.融合MODIS与ASTER数据生成高空间分辨率时间序列LAI方法研究.《北京师范大学学报(自然科学版)》.2007,第43卷(第3期),第303-308页. * |
| 地块支持下MODIS-NDVI时间序列冬小麦种植面积测量研究;李乐等;《光谱学与光谱分析》;20110531;第31卷(第5期);第l379-1383页 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN102831310A (en) | 2012-12-19 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN102831310B (en) | Method for building high-spatial resolution NDVI (normalized difference vegetation index) time series data | |
| Potapov et al. | Quantifying forest cover loss in Democratic Republic of the Congo, 2000–2010, with Landsat ETM+ data | |
| Molero et al. | SMOS disaggregated soil moisture product at 1 km resolution: Processor overview and first validation results | |
| Yang et al. | Urban land-cover change detection through sub-pixel imperviousness mapping using remotely sensed data | |
| Alonzo et al. | Mapping tall shrub biomass in Alaska at landscape scale using structure-from-motion photogrammetry and lidar | |
| Hubbard et al. | Comparative lahar hazard mapping at Volcan Citlaltépetl, Mexico using SRTM, ASTER and DTED-1 digital topographic data | |
| Fischer et al. | Monitoring topographic changes in a periglacial high‐mountain face using high‐resolution DTMs, Monte Rosa East Face, Italian Alps | |
| AU2012101249A4 (en) | Method for Generating High Spatial Resolution NDVI Time Series Data | |
| Conforti et al. | Morphological changes detection of a large earthflow using archived images, lidar-derived dtm, and uav-based remote sensing | |
| Hoffmann | Mapping damage during the Bam (Iran) earthquake using interferometric coherence | |
| Zandbergen | Applications of shuttle radar topography mission elevation data | |
| Berthier et al. | Measuring glacier mass changes from space—a review | |
| Liang et al. | A cross comparison of spatiotemporally enhanced springtime phenological measurements from satellites and ground in a northern US mixed forest | |
| Gillan et al. | Spatially explicit rangeland erosion monitoring using high-resolution digital aerial imagery | |
| Rogers et al. | Assessment of elevation uncertainty in salt marsh environments using discrete-return and full-waveform lidar | |
| Stângă et al. | Using old maps and soil properties to reconstruct the forest spatial pattern in the late 18th century. | |
| Strunk et al. | Evaluation of pushbroom DAP relative to frame camera DAP and lidar for forest modeling | |
| Ticehurst et al. | Using passive microwave and optical remote sensing to monitor flood inundation in support of hydrologic modelling | |
| Ziegler et al. | Using GEDI as training data for an ongoing mapping of landscape-scale dynamics of the plant area index | |
| Van Den Hoek et al. | Examining the utility of satellite-based wind sheltering estimates for lake hydrodynamic modeling | |
| Small et al. | Reliably flattened radar backscatter for wet snow mapping from wide-swath sensors | |
| Scarth et al. | Management effects on ground cover clumpiness: scaling from field to Sentinel-2 cover estimates | |
| Fry et al. | A low-cost GPS-based protocol to create high-resolution digital elevation models for remote mountain areas | |
| Bhushan et al. | Deriving seasonal and annual surface mass balance for debris-covered glaciers from flow-corrected satellite stereo DEM time series | |
| Jafari et al. | Comprehensive introduction to Digital Elevation Models, as a key dataset in soil erosion mapping |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| C14 | Grant of patent or utility model | ||
| GR01 | Patent grant | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20150401 Termination date: 20150817 |
|
| EXPY | Termination of patent right or utility model |