CN111768101B - Remote sensing cultivated land change detection method and system taking account of physical characteristics - Google Patents
Remote sensing cultivated land change detection method and system taking account of physical characteristics Download PDFInfo
- Publication number
- CN111768101B CN111768101B CN202010608516.XA CN202010608516A CN111768101B CN 111768101 B CN111768101 B CN 111768101B CN 202010608516 A CN202010608516 A CN 202010608516A CN 111768101 B CN111768101 B CN 111768101B
- Authority
- CN
- China
- Prior art keywords
- time sequence
- cultivated land
- model
- change
- data set
- 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.)
- Active
Links
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
- G06Q10/00—Administration; Management
- G06Q10/06—Resources, workflows, human or project management; Enterprise or organisation planning; Enterprise or organisation modelling
- G06Q10/063—Operations research, analysis or management
- G06Q10/0639—Performance analysis of employees; Performance analysis of enterprise or organisation operations
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
- G06Q10/00—Administration; Management
- G06Q10/06—Resources, workflows, human or project management; Enterprise or organisation planning; Enterprise or organisation modelling
- G06Q10/063—Operations research, analysis or management
- G06Q10/0637—Strategic management or analysis, e.g. setting a goal or target of an organisation; Planning actions based on goals; Analysis or evaluation of effectiveness of goals
- G06Q10/06375—Prediction of business process outcome or impact based on a proposed change
Abstract
The invention discloses a remote sensing cultivated land change detection method and system considering physical characteristics. The method comprises the following steps: constructing a plurality of time sequence data sets according to the time sequence vegetation indexes; processing each time sequence data set by adopting an improved MPFIT algorithm to obtain a processed time sequence data set; constructing a multi-harmonic model corresponding to each time sequence data set according to the processed time sequence data sets and the characteristics of double peaks or multiple peaks distribution in the physical years of the cultivated land; optimizing each multi-harmonic model to obtain a cultivated land feature time sequence track model corresponding to each time sequence data set; calculating the change intensity of the cultivated land based on the time sequence track model of each cultivated land object; judging whether the cultivated land changes according to the change intensity. The invention increases from original change detection based on two points to change detection based on two climatic time sequence tracks, improves the accuracy, precision and robustness of the change detection, and ensures the authenticity of the change information.
Description
Technical Field
The invention relates to the field of remote sensing application and analysis, in particular to a remote sensing cultivated land change detection method and system considering physical characteristics.
Background
Timely and accurately acquiring the spatial distribution of cultivated land and the change information thereof has great significance for optimizing and configuring cultivated land resources, formulating protection policies, researching ecological environment, planning sustainable development and the like. Remote sensing has been widely used for many years to determine spatial distribution and change analysis of large-scale cultivated land because of its macroscopic, efficient and convenient advantages. The remote sensing detection of the change of the cultivated land refers to the final determination and analysis of the change condition of the cultivated land, including the change position, the change range, the change type and other information, by utilizing remote sensing images covering the same area in different periods and combining with the change of the spectral characteristics of the cultivated land through a certain change detection method.
The traditional change detection method is mainly based on the spectrum difference of two images, such as a Change Vector Analysis (CVA), a post-classification comparison (PCC) and a correlation coefficient method, is simple to calculate and easy to understand, but is difficult to eliminate the ground feature pseudo-change phenomenon caused by random interference and different seasons. Especially in the detection of large-scale tilling changes, the false changes (false changes, i.e. all but tilling type changes) are more severe. The accumulation of massive historical remote sensing data enables more and more change detection methods based on time series images, but the existing research is mainly focused on forest disturbance, town expansion and other types, and a time series change detection method for cultivated lands is still lacking. The reason for this is mainly that due to the differences of crop types, growth stages, cultivation systems and the like, cultivated lands have more complicated heterogeneity in spectral characteristics and stronger dynamic properties in time sequence changes compared with forests and towns. Therefore, the conventional time-series change detection method is not effective for detecting the change of the cultivated land.
Disclosure of Invention
The invention aims to provide a remote sensing farmland change detection method and a remote sensing farmland change detection system taking account of physical characteristics, which are used for improving the accuracy, precision and robustness of change detection and ensuring the authenticity of change information from original change detection based on two points to change detection based on two physical time sequence tracks.
In order to achieve the above object, the present invention provides the following solutions:
a method for detecting changes in a remote sensing cultivated land in consideration of climatic features, the method comprising:
performing linear transformation on time sequence remote sensing image data acquired by a terrestrial satellite to obtain a time sequence vegetation index;
constructing a plurality of time sequence data sets according to the time sequence vegetation index, wherein the time sequence vegetation index is a first time sequence data set and a second time sequence data set;
processing each time sequence data set by adopting an improved MPFIT algorithm to obtain a processed time sequence data set;
constructing a multi-harmonic model corresponding to each time sequence data set according to each processed time sequence data set and the characteristics of double peaks or multiple peaks distribution in the cultivated land considered climate year;
optimizing each multi-harmonic model to obtain a cultivated land feature time sequence track model corresponding to each time sequence data set;
calculating the change intensity of the cultivated land based on each cultivated land physical time sequence track model;
judging whether the cultivated land changes according to the change intensity.
Optionally, the processing each time sequence data set by adopting the improved MPFIT algorithm to obtain a processed time sequence data set specifically includes:
and performing block and loop processing on each time sequence data set by adopting an improved MPFIT algorithm.
Optionally, optimizing each multi-harmonic model to obtain a cultivated land feature time sequence track model corresponding to each time sequence data set specifically includes:
calculating the error square sum, the decision coefficient and the root mean square error of each multi-harmonic model;
determining a maximum error value according to the error square sum, the decision coefficient and the root mean square error;
and removing the maximum error value to obtain a cultivated land physical time sequence track model corresponding to each time sequence data set.
Optionally, calculating the change intensity of the cultivated land based on each cultivated land feature time sequence track model specifically includes:
calculating coefficients of each cultivated land physical time sequence track model;
calculating the amplitude and the phase of each order of harmonic according to the coefficient;
and calculating the variation intensity according to the amplitude and the phase.
Optionally, after judging whether the cultivated land changes according to the change intensity, further including:
constructing a database by taking the reference classified images as a benchmark; the database comprises model coefficient ratio vectors;
calculating the minimum distance between the variation intensity and the model coefficient ratio vector;
and determining the change type according to the minimum distance.
The invention also provides a remote sensing cultivated land change detection system taking the physical characteristics into consideration, which comprises:
the linear transformation module is used for carrying out linear transformation on the time sequence remote sensing image data acquired by the terrestrial satellite to obtain a time sequence vegetation index;
the time sequence data set construction module is used for constructing a plurality of time sequence data sets which are a first time sequence data set and a second time sequence data set according to the time sequence vegetation index;
the data set processing module is used for processing each time sequence data set by adopting an improved MPFIT algorithm to obtain a processed time sequence data set;
the model construction module is used for constructing a multi-harmonic model corresponding to each time sequence data set according to each processed time sequence data set and the double-peak or double-peak distribution characteristics of the cultivated land in the object year;
the optimization module is used for optimizing each multi-harmonic model to obtain a cultivated land feature time sequence track model corresponding to each time sequence data set;
the change intensity calculation module is used for calculating the change intensity of the cultivated land based on each cultivated land feature time sequence track model;
and the cultivated land change judging module is used for judging whether the cultivated land changes according to the change intensity.
Optionally, the data set processing module adopts an improved MPFIT algorithm to perform block and loop processing on each time sequence data set.
Optionally, the optimizing module specifically includes:
the calculating unit is used for calculating the error square sum, the decision coefficient and the root mean square error of each multi-harmonic model;
a maximum error value determining unit, configured to determine a maximum error value according to the sum of squares of the errors, the decision coefficient, and the root mean square error;
and the removing unit is used for removing the maximum error value to obtain a cultivated land object time sequence track model corresponding to each time sequence data set.
Optionally, the variable intensity calculating module specifically includes:
the coefficient calculation unit is used for calculating the coefficient of each cultivated land feature time sequence track model;
an amplitude and phase calculation unit for calculating the amplitude and phase of each order harmonic according to the coefficient;
and the change intensity calculation unit is used for calculating the change intensity according to the amplitude and the phase.
Optionally, the method further comprises:
the data construction module is used for constructing a database by taking the reference classified images as the standard; the database comprises model coefficient ratio vectors;
the distance calculation module is used for calculating the minimum distance between the change intensity and the model coefficient ratio vector;
and the change type determining module is used for determining the change type according to the minimum distance.
According to the specific embodiment provided by the invention, the invention discloses the following technical effects:
(1) According to the invention, unique physical characteristics of cultivated land are considered, the change detection is increased from the original change detection based on two points to the change detection based on two physical time sequence tracks, the pseudo change caused by different seasons is eliminated, and the accuracy and the robustness of the extraction of the cultivated land change information are good.
(2) The invention can accurately detect real type change information, such as the change from cultivated land to built-up area, from cultivated land to water body, etc.; meanwhile, the false change caused by the natural state change of the cultivated land in the growing season can be well eliminated.
(3) The cultivated land change information extracted by the method can be used for analyzing the space-time change of the cultivated land, and further is used for evaluating grain yield, deciding resource allocation and the like.
Drawings
In order to more clearly illustrate the embodiments of the present invention or the technical solutions of the prior art, the drawings that are needed in the embodiments will be briefly described below, it being obvious that the drawings in the following description are only some embodiments of the present invention, and that other drawings may be obtained according to these drawings without inventive effort for a person skilled in the art.
FIG. 1 is a flow chart of a method for detecting changes in a remote sensing cultivated land according to an embodiment of the invention, which takes physical characteristics into consideration;
FIG. 2 is a schematic diagram of a multi-harmonic model construction taking account of the distribution of the weathered double (multi) peaks according to an embodiment of the present invention;
FIG. 3 is a schematic diagram of the calculation of the variation intensity and variation type according to the embodiment of the present invention;
fig. 4 is a block diagram of a remote sensing cultivated land change detection system according to an embodiment of the present invention, which takes physical characteristics into consideration.
Detailed Description
The following description of the embodiments of the present invention will be made clearly and completely with reference to the accompanying drawings, in which it is apparent that the embodiments described are only some embodiments of the present invention, but not all embodiments. All other embodiments, which can be made by those skilled in the art based on the embodiments of the invention without making any inventive effort, are intended to be within the scope of the invention.
The invention aims to provide a remote sensing farmland change detection method and a remote sensing farmland change detection system taking account of physical characteristics, which are used for improving the accuracy, precision and robustness of change detection and ensuring the authenticity of change information from original change detection based on two points to change detection based on two physical time sequence tracks.
In order that the above-recited objects, features and advantages of the present invention will become more readily apparent, a more particular description of the invention will be rendered by reference to the appended drawings and appended detailed description.
As shown in fig. 1, a remote sensing cultivated land change detection method considering the climatic characteristics includes:
step 101: and linearly transforming the time sequence remote sensing image data acquired by the terrestrial satellite to obtain a time sequence vegetation index.
Step 102: and constructing a plurality of time sequence data sets according to the time sequence vegetation index, wherein the time sequence data sets are a first time sequence data set and a second time sequence data set.
Adopts a pair ofVegetation index EVI with very sensitive cultivated land features is respectively constructed as t 1 And t 2 Two years worth of EVI timing feature data sets (i.e., a first timing data set and a second timing data set) are formulated as follows:
wherein R is blue ，R red ，R NIR Blue band, red band and near infrared band respectively representing Landsat data are all from the land satellite Landsat image data.
Step 103: and processing each time sequence data set by adopting an improved MPFIT algorithm to obtain a processed time sequence data set. Specifically, the improved MPFIT algorithm is adopted to carry out blocking and circulating processing on each time sequence data set.
While the MPFIT algorithm in IDL8.0 gives the best fit, it is currently only a fit for a single time series. Therefore, if the model coefficient estimation of each pixel on each time sequence index is performed according to the original MPFIT algorithm, a large memory space is required and the calculation amount is huge. Therefore, under the condition of ensuring the precision, in order to improve the calculation efficiency, the MPFIT algorithm is improved, so that the MPFIT algorithm can be subjected to blocking and circulating processing.
The modified MPFIT algorithm is processed in two steps, the first step: firstly, an ENVI_INIT_TILE function is adopted to block original data, then, after the function ENVI_GET_TILE is utilized to obtain block data, model fitting estimation is carried out through an MPFIT algorithm, and finally, the ENVI_TILE_DONE function is introduced to timely release the block data. And a second step of: on the basis of blocking processing, ENVI_BATCH_INIT is firstly introduced to set a BATCH processing inlet, then model coefficient estimation of each pixel is carried out by nesting 3 FOR … DOBEGIN loops, and finally BATCH loop estimation of model coefficients of each pixel on three time sequence indexes (NDVI, EVI and LSWI) is realized. Improving the MPFIT will find a set of coefficients that best fit the data under certain constraints. The blocking processing is to divide the original data of the whole scene EVI remote sensing image into smaller image blocks and then sequentially process the image blocks, so that the problem of insufficient memory in the big data processing process is solved. The cyclic processing is to automatically process the images or files of the same type in batches, so as to improve the calculation efficiency.
Step 104: and constructing a multi-harmonic model corresponding to each time sequence data set according to each processed time sequence data set and the characteristics of double-peak or multi-peak distribution in the physical year of the cultivated land, wherein the multi-harmonic model is shown in fig. 2.
Wherein t represents the acquisition days in the data year; t=365;
a represents the fitting value of the vegetation index at time t, a k Representing amplitude, b k Representing the phase.
Step 105: and optimizing each multi-harmonic model to obtain a cultivated land feature time sequence track model corresponding to each time sequence data set. Specific: calculating the error square sum, the decision coefficient and the root mean square error of each multi-harmonic model; determining a maximum error value according to the error square sum, the decision coefficient and the root mean square error; and removing the maximum error value to obtain a cultivated land physical time sequence track model corresponding to each time sequence data set.
The multi-harmonic model often involves the cloud and the pixel or noise point shielded by the cloud shadow into fitting, which results in poor model robustness. Therefore, in order to remove noise points that deviate from the normal range, a best fit model is preferable, and coefficients (R) are determined based on the Sum of Squares of Errors (SSE) 2 ) And Root Mean Square Error (RMSE) testing the stability and robustness of the model. When the model fitting accuracy is lower than a certain condition (R 2 < 0.6 and RMSE > 0.5 and SSE > 1.0), the observation point with the largest error value is removed and the remaining observations are re-fitted until a best fit model is obtained.
Step 106: and calculating the change intensity of the cultivated land based on each cultivated land physical time sequence track model. Specifically, calculating coefficients of each cultivated land physical time sequence track model; calculating the amplitude and the phase of each order of harmonic according to the coefficient; and calculating the variation intensity according to the amplitude and the phase.
And the variation intensity is calculated by adopting comprehensive similarity indexes based on amplitude, phase and residual error. The specific method comprises the following steps:
(1) After obtaining the weathered trajectory using the multiple harmonic model, 5 model coefficients (a 0 ，a 1 ，b 1 ，a 2 ，b 2 ) And the root mean square error RMSE of the model fitting are also obtained at the same time by extracting the model coefficient vector cv= (a) 0 ，a 1 ，b 1 ，a 2 ，b 2 ，RMSE) T And the climatic information hidden in the time sequence track is fully excavated, and the climatic track change trend is effectively expressed.
(2) Amplitude component A of each order harmonic k The expression is as follows:
(3) Phase component of each order harmonic
The expression is as follows: />
(4) Assume that
Respectively a certain pixel t 1 And t 2 Model coefficient vector for years, then the variation intensity CVD is calculated as follows:
step 107: judging whether the cultivated land changes according to the change intensity.
The difference image, namely the image with variable intensity is obtained through the steps. Theoretically, when the change intensity is 0, it is indicated that no change is found in the pixels for two years, and a value greater than 0 indicates that a change occurs. However, in practice, due to the influence of time phase differences, illumination, atmosphere and other factors in the remote sensing image, there is often a pseudo-change phenomenon, so that an optimal threshold value needs to be determined to represent a real change pixel, namely, extraction of change information. In order to accurately judge whether the cultivated land is changed or not, the influence of the pseudo change is removed, a change threshold delta is required to be determined, if CVD > delta, the pixel is judged to be a changed pixel, otherwise, the pixel is not changed. All the change pixels form a change area.
The determination of the change threshold relates to the accuracy of the judgment of the change of the cultivated land, and the embodiment of the invention adopts the maximum inter-class variance method proposed by Otsu and 1978 to determine the change threshold. The method is a global-based binarization algorithm, and divides an image into a foreground (changing) part and a background (unchanged) part according to the gray level characteristic of the image. The difference between the two parts should be the largest when the optimal threshold is taken, and the criterion used in the OTSU algorithm to measure the difference is the more common maximum inter-class variance. If the inter-class variance between the foreground and the background is larger, the difference between the two parts forming the image is larger, when a part of targets are divided into the background by mistake or a part of the background is divided into the targets by mistake, the difference between the two parts is smaller, and when the division of the taken threshold value maximizes the inter-class variance, the probability of the wrong division is minimum.
After step 107, further includes:
step 108: constructing a database by taking the reference classified images as a benchmark; the database comprises model coefficient ratio vectors;
step 109: calculating the minimum distance between the variation intensity and the model coefficient ratio vector;
step 1010: and determining the change type according to the minimum distance.
The change area on a scene image is extracted through step 106, but specific change types are not clear, such as detailed change information of whether the cultivated land becomes a water body or a built-up area. The identification of the type of change is based on the change region. In general, different ground classes have different climatic trajectories (EVI timing curves). Similarly, different types of changes may have specific rules or patterns of change. When one pixel changes from farmland to other types, it causes a change in CV. Based on this idea, the present invention proposes a method of model Coefficient Ratio Vector (CRV) for quantitatively and qualitatively determining the type of change. Firstly, constructing a CRV rule base of a reference change type by taking a reference classified image as a reference, then calculating CRVs of the change pixels, and finally determining the change type by utilizing the minimum distance between the CRVs of the reference change type and the CRVs of the change pixels. Experience statistics shows that the model coefficient vector keeps consistency in the changing direction and the magnitude of the model coefficient vector of the same ground object type, and in contrast, when the town and the water body are changed by cultivated land, the model coefficient vector is changed in different directions, and the magnitude of the model coefficient vector is obviously different. Thus, the present invention requires integrating the direction and magnitude of the vector to identify the type of change. To keep the vector direction, the gap between different types of variation is enlarged. The type of change is identified herein by the method of constructing Coefficient Ratio Vectors (CRV).
As shown in fig. 3, the specific method for identifying the change type of the present invention comprises:
let t be 1 Coefficient ratio vector between various ground object types on the image is t 1 ，t 2 The coefficient ratio vector corresponding to the period has equivalence. First at t 1 The time period image is a standard image, a certain number of samples are selected according to different ground object types, model coefficients of each ground object are extracted, then the average model coefficient of each ground object is calculated to be used as a reference model coefficient, and finally coefficient ratio vector CRV between different ground objects is calculated to be used as a reference
Wherein i and j respectively represent any two types of ground objects;
the average coefficient vectors of the i-ground class and the j-ground class are represented, respectively.
(2) The method of Coefficient Ratio Vector (CRV) is adopted to identify the change type. The CRV not only can eliminate the pseudo-change caused by crop type difference, climate period shift or advance, but also can well identify the change type of the change pixels. Taking the established reference change type CRV as a standard vector, and then calculating the pixel from t 1 To t 2 Coefficient ratio vector for epoch:
(3) Finally, assume that a pixel is at t 1 The period belongs to the cultivated land and the corresponding reference CRV from cultivated land to all other types is selected. Finally, calculate the change pixel
+.>
The distance D between them, the type of the change picture element is divided into the change type at the minimum distance.
As shown in fig. 4, the present invention further provides a remote sensing cultivated land change detection system considering the characteristics of the weather, the system comprising:
the linear transformation module 401 is configured to perform linear transformation on the time-sequence remote sensing image data acquired by the terrestrial satellite, so as to obtain a time-sequence vegetation index.
The time series data set constructing module 402 is configured to construct a plurality of time series data sets, which are a first time series data set and a second time series data set, according to the time series vegetation index.
The data set processing module 403 is configured to process each of the time-series data sets by using an improved MPFIT algorithm, so as to obtain a processed time-series data set. The data set processing module adopts an improved MPFIT algorithm to carry out block and circulation processing on each time sequence data set.
The model construction module 404 is configured to construct a multi-harmonic model corresponding to each time series data set according to each processed time series data set and the characteristics of dual-peak or dual-peak distribution in the cultivated land considered object year.
And the optimizing module 405 is configured to optimize each of the multiple harmonic models to obtain a cultivated land feature time sequence track model corresponding to each time sequence data set.
The optimizing module 405 specifically includes:
the calculating unit is used for calculating the error square sum, the decision coefficient and the root mean square error of each multi-harmonic model;
a maximum error value determining unit, configured to determine a maximum error value according to the sum of squares of the errors, the decision coefficient, and the root mean square error;
and the removing unit is used for removing the maximum error value to obtain a cultivated land object time sequence track model corresponding to each time sequence data set.
The change intensity calculation module 406 is configured to calculate the change intensity of the cultivated land based on each of the cultivated land feature time sequence track models.
The change intensity calculation module 406 specifically includes:
the coefficient calculation unit is used for calculating the coefficient of each cultivated land feature time sequence track model;
an amplitude and phase calculation unit for calculating the amplitude and phase of each order harmonic according to the coefficient;
and the change intensity calculation unit is used for calculating the change intensity according to the amplitude and the phase.
The cultivated land change judging module 407 is configured to judge whether the cultivated land changes according to the change strength.
The system further comprises:
a data construction module 408, configured to construct a database based on the reference classified images; the database includes model coefficient ratio vectors.
A distance calculating module 409, configured to calculate a minimum distance between the variation intensity and the model coefficient ratio vector.
The change type determining module 4010 is configured to determine a change type according to the minimum distance.
In the present specification, each embodiment is described in a progressive manner, and each embodiment is mainly described in a different point from other embodiments, and identical and similar parts between the embodiments are all enough to refer to each other. For the system disclosed in the embodiment, since it corresponds to the method disclosed in the embodiment, the description is relatively simple, and the relevant points refer to the description of the method section.
The principles and embodiments of the present invention have been described herein with reference to specific examples, the description of which is intended only to assist in understanding the methods of the present invention and the core ideas thereof; also, it is within the scope of the present invention to be modified by those of ordinary skill in the art in light of the present teachings. In view of the foregoing, this description should not be construed as limiting the invention.
Claims (8)
1. A method for detecting changes in a remote sensing cultivated land in consideration of physical characteristics, the method comprising:
performing linear transformation on time sequence remote sensing image data acquired by a terrestrial satellite to obtain a time sequence vegetation index;
constructing a plurality of time sequence data sets according to the time sequence vegetation index, wherein the time sequence vegetation index is a first time sequence data set and a second time sequence data set;
processing each time sequence data set by adopting an improved MPFIT algorithm to obtain a processed time sequence data set; the method specifically comprises the following steps: performing block and circulation processing on each time sequence data set by adopting an improved MPFIT algorithm; the modified MPFIT algorithm is processed in two steps, the first step: firstly, an ENVI_INIT_TILE function is adopted to block original data, then, after the block data are obtained by utilizing the ENVI_GET_TILE function, model fitting estimation is carried out through an MPFIT algorithm, and finally, the ENVI_TILE_DONE function is introduced to timely release the block data; and a second step of: on the basis of blocking processing, firstly introducing ENVI_BATCH_INIT to set a BATCH processing inlet, then carrying out model coefficient estimation of each pixel by nesting 3 FOR … DO BEGIN loops, and finally realizing BATCH loop estimation of model coefficients of each pixel on three time sequence indexes NDVI, EVI and LSWI; wherein, NDVI is normalized vegetation index, EVI is enhanced vegetation index, LSWI is vegetation moisture content index;
constructing a multi-harmonic model corresponding to each time sequence data set according to each processed time sequence data set and the characteristics of double peaks or multiple peaks distribution in the cultivated land considered climate year; the expression of the multi-harmonic model is as follows:
wherein t represents the acquisition days in the data year; t=365;
a represents the fitting value of the vegetation index at time t, a k Representing amplitude, b k Representing phase;
optimizing each multi-harmonic model to obtain a cultivated land feature time sequence track model corresponding to each time sequence data set;
calculating the change intensity of the cultivated land based on each cultivated land physical time sequence track model; the variation intensity CVD was calculated as follows:
wherein a is k Representing amplitude, b k The phase is represented by a phase value,
、/>
respectively a certain pixel t 1 And t 2 The model coefficient vector for the year,
、/>
respectively a certain pixel t 1 And t 2 Annual model root mean square error; />
、/>
Respectively a certain pixel t 1 And t 2 Amplitude of each order harmonic of the annual multiple harmonic model; />
、/>
Respectively a certain pixel t 1 And t 2 The phase of each order of harmonics of the annual multiple harmonic model;
judging whether the cultivated land changes according to the change intensity.
2. The method for detecting the change of the cultivated land by remote sensing according to claim 1, wherein the optimizing each multi-harmonic model to obtain a cultivated land physical time sequence track model corresponding to each time sequence data set specifically comprises:
calculating the error square sum, the decision coefficient and the root mean square error of each multi-harmonic model;
determining a maximum error value according to the error square sum, the decision coefficient and the root mean square error;
and removing the maximum error value to obtain a cultivated land physical time sequence track model corresponding to each time sequence data set.
3. The method for detecting a change in a remote sensing cultivated land according to claim 1, wherein the calculating the change intensity of the cultivated land based on each of the cultivated land feature time series track models comprises:
calculating coefficients of each cultivated land physical time sequence track model;
calculating the amplitude and the phase of each order of harmonic according to the coefficient;
calculating a variation intensity from the amplitude and the phase;
wherein the amplitude A of each order of harmonic k The expression is as follows:
；
the phase of each order of harmonic
The expression is as follows: />
。
4. The method for detecting a change in a remote sensing cultivated land in consideration of a physical characteristic according to claim 1, further comprising, after determining whether the cultivated land has changed according to the change intensity:
constructing a database by taking the reference classified images as a benchmark; the database comprises model coefficient ratio vectors;
calculating the minimum distance between the variation intensity and the model coefficient ratio vector;
determining a change type according to the minimum distance;
the method specifically comprises the following steps:
calculating coefficient ratio vector CRV between different ground objects as reference
；
Wherein i and j respectively represent any two types of ground objects;
，/>
average coefficient vectors respectively representing i-ground class and j-ground class; />
，/>
，/>
，/>
，/>
5 model coefficient vectors respectively representing each pixel of the i-field class; />
，/>
，/>
，/>
，/>
5 model coefficient vectors respectively representing each pixel of j-th class；/>
，/>
Model root mean square error vectors respectively representing i-place class and j-place class;
identifying the change type by adopting a coefficient ratio vector CRV method, taking the established reference change type CRV as a standard vector, and then calculating the pixel from t 1 To t 2 Coefficient ratio vector for epoch:
calculating a change pixel
+.>
Distance D between, the type of the change pixel is divided into the change type at the minimum distance:
5. a remote sensing cultivated land change detection system taking into account climatic features, the system comprising:
the linear transformation module is used for carrying out linear transformation on the time sequence remote sensing image data acquired by the terrestrial satellite to obtain a time sequence vegetation index;
the time sequence data set construction module is used for constructing a plurality of time sequence data sets which are a first time sequence data set and a second time sequence data set according to the time sequence vegetation index;
the data set processing module is used for processing each time sequence data set by adopting an improved MPFIT algorithm to obtain a processed time sequence data set; the data set processing module adopts an improved MPFIT algorithm to carry out block and circulation processing on each time sequence data set; the modified MPFIT algorithm is processed in two steps, the first step: firstly, an ENVI_INIT_TILE function is adopted to block original data, then, after the block data are obtained by utilizing the ENVI_GET_TILE function, model fitting estimation is carried out through an MPFIT algorithm, and finally, the ENVI_TILE_DONE function is introduced to timely release the block data; and a second step of: on the basis of blocking processing, firstly introducing ENVI_BATCH_INIT to set a BATCH processing inlet, then carrying out model coefficient estimation of each pixel by nesting 3 FOR … DO BEGIN loops, and finally realizing BATCH loop estimation of model coefficients of each pixel on three time sequence indexes NDVI, EVI and LSWI; wherein, NDVI is normalized vegetation index, EVI is enhanced vegetation index, LSWI is vegetation moisture content index;
the model construction module is used for constructing a multi-harmonic model corresponding to each time sequence data set according to each processed time sequence data set and the double-peak or double-peak distribution characteristics of the cultivated land in the object year; the expression of the multi-harmonic model is as follows:
wherein t represents the acquisition days in the data year; t=365;
a represents the fitting value of the vegetation index at time t, a k Representing amplitude, b k Representing phase;
the optimization module is used for optimizing each multi-harmonic model to obtain a cultivated land feature time sequence track model corresponding to each time sequence data set;
the change intensity calculation module is used for calculating the change intensity of the cultivated land based on each cultivated land feature time sequence track model; the variation intensity CVD was calculated as follows:
wherein a is k Representing amplitude, b k The phase is represented by a phase value,
、/>
respectively a certain pixel t 1 And t 2 The model coefficient vector for the year,
、/>
respectively a certain pixel t 1 And t 2 Annual model root mean square error; />
、/>
Respectively a certain pixel t 1 And t 2 Amplitude of each order harmonic of the annual multiple harmonic model; />
、/>
Respectively a certain pixel t 1 And t 2 The phase of each order of harmonics of the annual multiple harmonic model;
and the cultivated land change judging module is used for judging whether the cultivated land changes according to the change intensity.
6. The system for detecting changes in a remote sensing cultivated land taking account of climatic features according to claim 5, wherein said optimization module comprises:
the calculating unit is used for calculating the error square sum, the decision coefficient and the root mean square error of each multi-harmonic model;
a maximum error value determining unit, configured to determine a maximum error value according to the sum of squares of the errors, the decision coefficient, and the root mean square error;
and the removing unit is used for removing the maximum error value to obtain a cultivated land object time sequence track model corresponding to each time sequence data set.
7. The system for detecting changes in a remote sensing cultivated land taking account of climatic features according to claim 5, wherein the change intensity calculation module specifically comprises:
the coefficient calculation unit is used for calculating the coefficient of each cultivated land feature time sequence track model;
an amplitude and phase calculation unit for calculating the amplitude and phase of each order harmonic according to the coefficient;
a change intensity calculation unit configured to calculate a change intensity from the amplitude and the phase;
wherein the amplitude A of each order of harmonic k The expression is as follows:
；/>
the phase of each order of harmonic
The expression is as follows: />
。
8. The system for detecting changes in a remote sensing cultivated land in consideration of physical characteristics according to claim 5, further comprising:
the data construction module is used for constructing a database by taking the reference classified images as the standard; the database comprises model coefficient ratio vectors;
the distance calculation module is used for calculating the minimum distance between the change intensity and the model coefficient ratio vector;
the change type determining module is used for determining a change type according to the minimum distance;
the method specifically comprises the following steps:
calculating coefficient ratio vector CRV between different ground objects as reference
；
Wherein i and j respectively represent any two types of ground objects;
，/>
average coefficient vectors respectively representing i-ground class and j-ground class; />
，/>
，/>
，/>
，/>
5 model coefficient vectors respectively representing each pixel of the i-field class; />
，/>
，/>
，/>
，/>
5 model coefficient vectors respectively representing each pixel of the j-land class; />
，/>
Model root mean square error vectors respectively representing i-place class and j-place class;
identifying the change type by adopting a coefficient ratio vector CRV method, taking the established reference change type CRV as a standard vector, and then calculating the pixel from t 1 To t 2 Coefficient ratio vector for epoch:
calculating a change pixel
+.>
Distance D between, the type of the change pixel is divided into the change type at the minimum distance:
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010608516.XA CN111768101B (en) | 2020-06-29 | 2020-06-29 | Remote sensing cultivated land change detection method and system taking account of physical characteristics |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010608516.XA CN111768101B (en) | 2020-06-29 | 2020-06-29 | Remote sensing cultivated land change detection method and system taking account of physical characteristics |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN111768101A CN111768101A (en) | 2020-10-13 |
| CN111768101B true CN111768101B (en) | 2023-05-23 |
Family
ID=72724112
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202010608516.XA Active CN111768101B (en) | 2020-06-29 | 2020-06-29 | Remote sensing cultivated land change detection method and system taking account of physical characteristics |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN111768101B (en) |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN113538388B (en) * | 2021-07-23 | 2022-10-11 | 中国电子科技集团公司第五十四研究所 | Arable land loss assessment method based on MODIS NDVI time sequence data |
Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103700092A (en) * | 2013-12-04 | 2014-04-02 | 中国科学院遥感与数字地球研究所 | Forest burned area automatic extraction method based on time sequence remote sensing image |
| CN104318270A (en) * | 2014-11-21 | 2015-01-28 | 东北林业大学 | Land cover classification method based on MODIS time series data |
| CN106021653A (en) * | 2016-05-06 | 2016-10-12 | 华南农业大学 | An NDVI time sequence reconstruction method and system |
| CN106355143A (en) * | 2016-08-25 | 2017-01-25 | 中国农业大学 | Seed maize field identification method and system based on multi-source and multi-temporal high resolution remote sensing data |
| CN107403157A (en) * | 2017-07-28 | 2017-11-28 | 中国科学院东北地理与农业生态研究所 | Region large scale crops pattern of farming extracting method based on MODIS data |
| CN109115770A (en) * | 2018-06-14 | 2019-01-01 | 中科禾信遥感科技（苏州）有限公司 | A kind of a wide range of crops remote-sensing monitoring method and device |
| CN110766299A (en) * | 2019-10-11 | 2020-02-07 | 清华大学 | Watershed vegetation change analysis method based on remote sensing data |
-
2020
- 2020-06-29 CN CN202010608516.XA patent/CN111768101B/en active Active
Patent Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103700092A (en) * | 2013-12-04 | 2014-04-02 | 中国科学院遥感与数字地球研究所 | Forest burned area automatic extraction method based on time sequence remote sensing image |
| CN104318270A (en) * | 2014-11-21 | 2015-01-28 | 东北林业大学 | Land cover classification method based on MODIS time series data |
| CN106021653A (en) * | 2016-05-06 | 2016-10-12 | 华南农业大学 | An NDVI time sequence reconstruction method and system |
| CN106355143A (en) * | 2016-08-25 | 2017-01-25 | 中国农业大学 | Seed maize field identification method and system based on multi-source and multi-temporal high resolution remote sensing data |
| CN107403157A (en) * | 2017-07-28 | 2017-11-28 | 中国科学院东北地理与农业生态研究所 | Region large scale crops pattern of farming extracting method based on MODIS data |
| CN109115770A (en) * | 2018-06-14 | 2019-01-01 | 中科禾信遥感科技（苏州）有限公司 | A kind of a wide range of crops remote-sensing monitoring method and device |
| CN110766299A (en) * | 2019-10-11 | 2020-02-07 | 清华大学 | Watershed vegetation change analysis method based on remote sensing data |
Non-Patent Citations (13)
| Title |
|---|
| A rule-based approach for crop identification using multi-temporaland multi-sensor phenological metrics;Taylor & Francia;第51卷(第1期);511–524 * |
| Agriculture Phenology Monitoring Using NDVI Time Series Based on Remote Sensing Satellites: A Case Study of Guangdong, China;Komal Choudhary, Wenzhong Shi, Mukesh Singh Boori & Samuel Corgne;Optical memory & neural networks;第28卷(第3期);204-214 * |
| Detection of Cropland Change Using Multi-Harmonic Based Phenological Trajectory Similarity;Jiage Chen ,Jun Chen ,Huiping Liu , andShu Peng;remoting sensing;第10卷;1-20 * |
| Determination of Vegetation Thresholds for Assessing Land Use and Land Use Changes in Cambodia using the Google Earth Engine Cloud-Computing Platform .remote sensing .2019,第11卷1-30. * |
| Land–Use and Land-Cover Change Detection Using Dynamic Time Warping–Based Time Series Clustering Method;Yanghua Zhang &Hu Zhao;Canadian Journal of Remote Sensing;第46卷(第1期);67-83 * |
| Long-Term Satellite Image Time-Series for Land Use/Land Cover Change Detection Using Refined Open Source Data in a Rural Region;Cláudia M. Viana *ORCID,Inês GirãoORCID and Jorge Rocha;Remote Sensing;第11卷(第9期);1-22 * |
| Spatio-Temporal Reconstruction of MODIS NDVI by Regional Land Surface Phenology and Harmonic Analysis of Time-Series;remoting sensing;第56卷(第8期);1261-1288 * |
| 农业土地利用遥感信息提取的研究进展与展;董金玮;吴文斌;黄健熙;尤南山;何盈利;地球信息科学学报;第22卷(第4期);772-783 * |
| 基于HJ卫星数据与面向对象分类的土地利用/覆盖信息提取;朱永森，曾永年,张猛;农业工程学报;第33卷(第14期);258-265 * |
| 基于时序植被指数的湖北省物候空间特征分析;李彩霞,邓帆,张佳华等;长江流域资源与环境;第28卷(第7期);1583-1589 * |
| 基于遥感时序分析的黑龙江省撂荒地时空监测;邓欣雨;中国优秀硕士学位论文全文数据库 (基础科学辑)(第2期);1-59 * |
| 李彩霞,邓帆,张佳华等.基于时序植被指数的湖北省物候空间特征分析.长江流域资源与环境.2019,第28卷(第7期),1583-1589. * |
| 线性内插和扩展Kalman滤波组合法对NDVI时间序列的重构与预测;蒋雪冰;中国优秀硕士学位论文全文数据库 (基础科学辑)(第3期);1-53 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN111768101A (en) | 2020-10-13 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN108921885B (en) | Method for jointly inverting forest aboveground biomass by integrating three types of data sources | |
| CN108846832B (en) | Multi-temporal remote sensing image and GIS data based change detection method and system | |
| CN106780091B (en) | Agricultural disaster information remote sensing extraction method based on vegetation index time-space statistical characteristics | |
| Wei et al. | Cloud detection for Landsat imagery by combining the random forest and superpixels extracted via energy-driven sampling segmentation approaches | |
| CN111598045B (en) | Remote sensing farmland change detection method based on object spectrum and mixed spectrum | |
| CN104751478B (en) | Object-oriented building change detection method based on multi-feature fusion | |
| CN102254319B (en) | Method for carrying out change detection on multi-level segmented remote sensing image | |
| CN103632363B (en) | Object level high-resolution remote sensing image change detecting method based on Multiscale Fusion | |
| CN112164062A (en) | Wasteland information extraction method and device based on remote sensing time sequence analysis | |
| CN113221765B (en) | Vegetation phenological period extraction method based on digital camera image effective pixels | |
| CN107330875B (en) | Water body surrounding environment change detection method based on forward and reverse heterogeneity of remote sensing image | |
| CN112381013B (en) | Urban vegetation inversion method and system based on high-resolution remote sensing image | |
| CN105488770A (en) | Object-oriented airborne laser radar point cloud filtering method | |
| CN108710864B (en) | Winter wheat remote sensing extraction method based on multi-dimensional identification and image noise reduction processing | |
| Lavanya et al. | Terrain mapping of LandSat8 images using MNF and classifying soil properties using ensemble modelling | |
| CN110765962A (en) | Plant identification and classification method based on three-dimensional point cloud contour dimension values | |
| Hu et al. | Integrating CART algorithm and multi-source remote sensing data to estimate sub-pixel impervious surface coverage: a case study from Beijing Municipality, China | |
| CN111091079A (en) | TLS-based method for measuring dominant single plant structural parameters of vegetation in alpine and fragile regions | |
| CN111339948A (en) | Automatic identification method for newly-added buildings of high-resolution remote sensing images | |
| CN115468917A (en) | Method and system for extracting crop information of farmland plot based on high-resolution remote sensing | |
| CN111768101B (en) | Remote sensing cultivated land change detection method and system taking account of physical characteristics | |
| He et al. | ICESat-2 data classification and estimation of terrain height and canopy height | |
| CN102231190B (en) | Automatic extraction method for alluvial-proluvial fan information | |
| CN117075138A (en) | Remote sensing measurement and calculation method, system and medium for canopy height of 30-meter forest in area | |
| Moradi et al. | Potential evaluation of visible-thermal UAV image fusion for individual tree detection based on convolutional neural network |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |