CN104616265A - Time domain reconstruction method of remote sensing sequence data - Google Patents

Abstract

The invention discloses a time domain reconstruction method of remote sensing sequence data. The time domain reconstruction method comprises the following steps: firstly, calculating the correlation of non-missing public parts between the remote sensing sequence data, and reconstructing reasonable distribution weights for missed parts through the correlation size; secondly, calculating a matching corresponding value of a reference image and an image to be reconstructed through moment matching; and finally, reconstructing a time sequence by a weighting harmonic analyzing method. By the aid of the time domain reconstruction method, a gray level of the reference image horizontally gets closer to the gray level of the image to be reconstructed by a moment matching principle through the correlation reasonable distribution weights between the data; time domain reconstruction of the remote sensing sequence data is carried out by the weighting harmonic analyzing method through abiding by a time dependence relationship, so that the time dependence relationship of the time sequence data is sufficiently utilized and a reconstructed result is real and available; the time domain reconstruction method is easy to realize and has high efficiency and important actual application meanings.

Description

A kind of time domain method for reconstructing of remote sensing sequence data

Technical field

The invention belongs to remote sensing image processing technology field, relate to a kind of time domain method for reconstructing of remote sensing sequence data, be specifically related to the remote sensing time series data time domain method for reconstructing of a kind of associating weighting frequency analysis and match by moment.

Background technology

Remote sensing long-term sequence data are widely used in global environmental change research, but pollute and many-sided impact such as uncertainty of air owing to being subject to cloud, cause remote sensing image and may there is a large amount of noise or loss of data, the normal use of remote sensing image is hindered.Therefore, develop rational method for reconstructing, for remote sensing time series data collection carries out denoising, repair process after acquisition, the quality and the stability that improve data are far reaching problems.

At present, the reparation for remote sensing image mainly divides based on single width image with based on several images two kinds.Based on the reparation of single width image, mainly rely on the information of self from mending the angle of painting, its main method has the method for interpolation, Histogram Matching, partial differential equation and total variation method etc.When repairing single width image, if disappearance area is less, can be repaired result preferably, during if there is large area information dropout, the information of single width image self is not enough to support the reconstruction to lack part, is difficult to well be repaired result.Based on the reparation of several images, compensate for the phenomenon of spatial information deficiency, conventional method has the method such as Histogram Matching and regretional analysis, data reconstruction can be carried out by setting up funtcional relationship between different images, if different images spectral characteristic exists very large difference, it is repaired result and often occurs significantly rebuilding vestige.At present, the process for long-term sequence data carries out denoising by the mode of filtering to data, and this little filtering method does not often consider the utilization of spatial information, and the data for consecutive miss are difficult to the reconstructed results obtained.

Summary of the invention

The present invention is directed to the shortcoming of prior art to space-time complementary information underutilization in time series, the remote sensing time series data time domain method for reconstructing of a kind of associating weighting frequency analysis and match by moment is proposed, on the basis utilizing the long-term sequence data time relations of dependence, consider spatial coherence, making full use of data space information in conjunction with match by moment theory, by adopting weighting harmonic analysis method, time domain reconstruction can be carried out to missing data accurately, and counting yield is high, practical.

The technical solution adopted in the present invention is: a kind of time domain method for reconstructing of remote sensing sequence data, is characterized in that, comprise the following steps:

Step 1: remote sensing time series data is carried out geometrical registration, obtains accurate registration image;

Step 2: carry out cloud detection to remote sensing time series data, obtains the cloud mask of each time data, and calculates the related coefficient CC of the public cloud-free area of each time series data;

Step 3: with reference to data public cloudless region grey level's approximate transform to the public cloudless region grey level of data to be reconstructed, generate new time series data;

Step 4: the related coefficient according to step 2 gained distributes weights, rebuilds new time series data.

As preferably, when step 2 carries out cloud detection, be also marked as cloud sector by the exceptional value of threshold decision.

As preferably, the related coefficient CC of the public cloud-free area of each time series data of the calculating described in step 2, its computing formula is:

$\mathrm{CC}(x,y)=\frac{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(x_i-m_x)(y_i-m_y)}{\sqrt{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{(x_i-m_x)}^2}\sqrt{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{(y_i-m_y)}^2}};$

Wherein x, y are two remotely-sensed datas, and N is the public cloud-free area number of pixels of remotely-sensed data, x _i, y _ii-th pixel in data x, y, 1≤i≤N, m _x, m _ythe average of the corresponding public cloud-free area of remotely-sensed data respectively.

As preferably, the specific implementation process of step 3 is, by moment-matching method with reference to the horizontal approximate transform of data gray value to data gray value level to be reconstructed, calculate a new time proximity sequence data in each data moment to be reconstructed; For moment t to be reconstructed _othe gray average of data is m _to, gray variance is s _to, other moment t _fthe gray average of data is m _tf, gray variance is s _tf, gray-scale value is g _tf, then new time proximity sequence data corresponding grey scale value is:

$g_{\mathrm{new}}=\frac{[g_{\mathrm{tf}}-m_{\mathrm{tf}}]s_{\mathrm{to}}}{s_{\mathrm{tf}}}+m_{\mathrm{to}}.$

As preferably, the specific implementation of step 4 comprises following sub-step:

Step 4.1: to reconstruction moment t _odata are recovered;

First be the new time series g of step 3 gained _newin other time data distribute weights according to the related coefficient (CC) of step 2 gained, its weight computing mode is:

$w_i=\frac{{\mathrm{CC}}_i}{\underset{i=1}{\overset{q}{\mathrm{\Σ}}}{\mathrm{CC}}_i};$

Wherein, q is the number without missing data in time series, CC _iit is the related coefficient of i-th data and the public cloud-free area of data to be reconstructed;

Next converts data to the form that sine and cosine is added, and answers gray-scale value g to the new time proximity sequence pair of step 3 gained _newbe weighted and solve, its implementation is:

$\begin{array}{c}\left|\right|g\left(t\right)-[a_0+\underset{j=1}{\overset{\mathrm{nf}}{\mathrm{\Σ}}}[a_j\cos \left(2\mathrm{\π}{f}_{j}t\right)+b_j\sin \left(2\mathrm{\π}{f}_{j}t\right)\left]\right]{\left|\right|}_w\\ =\∫[g\left(t_i\right)-[a_0+\underset{j=1}{\overset{\mathrm{nf}}{\mathrm{\Σ}}}{[a_0+\underset{j=1}{\overset{\mathrm{nf}}{\mathrm{\Σ}}}[a_j\cos \left(2\mathrm{\π}{f}_{j}t\right)+b_j\sin \left(2\mathrm{\π}{f}_{j}t\right)\left]\right]]}^2{w}_i\mathrm{dy}\end{array},;$

Wherein, nf is harmonic wave number, a ₀for coefficient when harmonic wave number is 0, a _jand b _jfor the coefficient of trigonometric function each several part, f _jfor frequency, t _jfor timing node, y is time proximity sequence data;

Step 4.2: travel through all single picture point time sequences, repeating above-mentioned steps can rebuild remote sensing time series data, the reconstruction mark problem caused with the grey level's difference eliminated between different images.

The remote sensing time series data time domain method for reconstructing of the associating weighting frequency analysis that the present invention proposes and match by moment utilizes the message complementary sense of remote sensing long-term sequence data to carry out time domain reconstruction.Not only take into account the temporal relations of dependence, and by the correlativity reasonable distribution weights between image, convert image greyscale level by moment-matching method, construct new time series and missing values is repaired, the reconstruction mark problem caused with the grey level's difference eliminated between different images.In a word, the method that the present invention proposes effectively can eliminate the impact of the factors such as cloud, improves the q&r of remote sensing time series data, for subsequent applications provides reliable support.Therefore, the method that the present invention proposes has important practical application meaning.

Accompanying drawing explanation

Fig. 1 is the process flow diagram of the embodiment of the present invention.

Embodiment

Understand for the ease of those of ordinary skill in the art and implement the present invention, below in conjunction with drawings and Examples, the present invention is described in further detail, should be appreciated that exemplifying embodiment described herein is only for instruction and explanation of the present invention, is not intended to limit the present invention.

Ask for an interview Fig. 1, the time domain method for reconstructing of a kind of remote sensing sequence data of the embodiment of the present invention, comprises the following steps:

Step 1: first need to carry out geometry correction to remote sensing time series data, obtain accurate registration image, concrete geometry correction is prior art.

Step 2: carry out cloud detection to remote sensing time series data, obtains the cloud mask of each time data, and calculates the related coefficient CC of the public cloud-free area of each time series data.

If remote sensing time series data exists quality status stamp wave band, carry out cloud Covering judgment by this wave band, if there is no this wave band then by other data or cloud detection algorithm carry out cloud cover detection set up cloud mask.Meanwhile, the exceptional value existed in data or invalid value are also labeled as cloud, are rebuild in subsequent treatment.After cloud detection, the related coefficient of sequence data computing time public cloud-free area each other, in order to reasonable distribution weights.Wherein the related coefficient of the public cloud-free area of each time series data is calculated by following formula:

$\mathrm{CC}(x,y)=\frac{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(x_i-m_x)(y_i-m_y)}{\sqrt{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{(x_i-m_x)}^2}\sqrt{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{(y_i-m_y)}^2}};$

Wherein x, y are two remotely-sensed datas, and N is the public cloud-free area number of pixels of remotely-sensed data, x _i, y _ii-th pixel in data x, y, 1≤i≤N, m _x, m _ythe average of the corresponding public cloud-free area of remotely-sensed data respectively.

Step 3: with reference to data public cloudless region grey level's approximate transform to the public cloudless region grey level of data to be reconstructed, generate new time series data.

It is process by picture point time series processing that this step is rebuild remote sensing time series data.Needed to carry out greyscale transformation in this pixel without missing data before execution is by picture point time series processing, generate a new time series.This process carries out match by moment by single picture point time sequence without missing data place image and this time series data place to be reconstructed image, and by the grey level of the horizontal approximate transform of nothing disappearance data gray to data place to be reconstructed image, its settling mode is:

$g_{\mathrm{new}}=\frac{[g_{\mathrm{tf}}-m_{\mathrm{tf}}]s_{\mathrm{to}}}{s_{\mathrm{tf}}}+m_{\mathrm{to}};$

Wherein, moment t to be reconstructed _othe cloud-free area gray average of data is m _to, gray variance is s _to, other moment t _fdata cloud-free area gray average is m _tf, gray variance is s _tf, gray-scale value is g _tf.

Step 4: distribute weights according to related coefficient, new time series is rebuild.

This step is carried out by picture point time series processing to remote sensing time series data, needs to rebuild missing data in new time series, to reconstruction moment t _odata are recovered, and utilize the method for weighting frequency analysis, convert data to the form that sine and cosine is added; Specific implementation comprises following sub-step:

Step 4.1: to reconstruction moment t _odata are recovered;

First be the new time series g of step 3 gained _newin other time data distribute weights according to the related coefficient (CC) of step 2 gained, its weight computing mode is:

$w_i=\frac{{\mathrm{CC}}_i}{\underset{i=1}{\overset{q}{\mathrm{\Σ}}}{\mathrm{CC}}_i};$

Wherein, q is the number without missing data in time series, CC _iit is the related coefficient of i-th data and the public cloud-free area of data to be reconstructed;

Next converts data to the form that sine and cosine is added, and answers gray-scale value g to the new time proximity sequence pair of step 3 gained _newbe weighted and solve, its implementation is:

$\begin{array}{c}\left|\right|g\left(t\right)-[a_0+\underset{j=1}{\overset{\mathrm{nf}}{\mathrm{\Σ}}}[a_j\cos \left(2\mathrm{\π}{f}_{j}t\right)+b_j\sin \left(2\mathrm{\π}{f}_{j}t\right)\left]\right]{\left|\right|}_w\\ =\∫[g\left(t_i\right)-[a_0+\underset{j=1}{\overset{\mathrm{nf}}{\mathrm{\Σ}}}{[a_0+\underset{j=1}{\overset{\mathrm{nf}}{\mathrm{\Σ}}}[a_j\cos \left(2\mathrm{\π}{f}_{j}t\right)+b_j\sin \left(2\mathrm{\π}{f}_{j}t\right)\left]\right]]}^2{w}_i\mathrm{dy}\end{array},;$

Wherein, nf is harmonic wave number, a ₀for coefficient when harmonic wave number is 0, a _jand b _jfor the coefficient of trigonometric function each several part, f _jfor frequency, t _jfor timing node, y is time proximity sequence data;

Step 4.2: travel through all single picture point time sequences, repeating above-mentioned steps can rebuild remote sensing time series data, the reconstruction mark problem caused with the grey level's difference eliminated between different images.

Should be understood that, the part that this instructions does not elaborate all belongs to prior art.

Should be understood that; the above-mentioned description for preferred embodiment is comparatively detailed; therefore the restriction to scope of patent protection of the present invention can not be thought; those of ordinary skill in the art is under enlightenment of the present invention; do not departing under the ambit that the claims in the present invention protect; can also make and replacing or distortion, all fall within protection scope of the present invention, request protection domain of the present invention should be as the criterion with claims.

Claims (5)

1. a time domain method for reconstructing for remote sensing sequence data, is characterized in that, comprise the following steps:

Step 1: remote sensing time series data is carried out geometrical registration, obtains accurate registration image;

Step 2: carry out cloud detection to remote sensing time series data, obtains the cloud mask of each time data, and calculates the related coefficient CC of the public cloud-free area of each time series data;

Step 3: with reference to data public cloudless region grey level's approximate transform to the public cloudless region grey level of data to be reconstructed, generate new time series data;

Step 4: the related coefficient according to step 2 gained distributes weights, rebuilds new time series data.

2. the time domain method for reconstructing of remote sensing sequence data according to claim 1, is characterized in that: when step 2 carries out cloud detection, is also marked as cloud sector by the exceptional value of threshold decision.

3. the time domain method for reconstructing of remote sensing sequence data according to claim 1, is characterized in that: the related coefficient CC of the public cloud-free area of each time series data of the calculating described in step 2, and its computing formula is:

$\mathrm{CC}(x,y)=\frac{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(x_i-m_x)(y_i-m_y)}{\sqrt{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{(x_i-m_x)}^2}\sqrt{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{(y_i-m_y)}^2}};$

Wherein x, y are two remotely-sensed datas, and N is the public cloud-free area number of pixels of remotely-sensed data, x _i, y _ii-th pixel in data x, y, 1≤i≤N, m _x, m _ythe average of the corresponding public cloud-free area of remotely-sensed data respectively.

4. the time domain method for reconstructing of remote sensing sequence data according to claim 3, it is characterized in that: the specific implementation process of step 3 is, by moment-matching method with reference to the horizontal approximate transform of data gray value to data gray value level to be reconstructed, calculate a new time proximity sequence data in each data moment to be reconstructed; For moment t to be reconstructed _othe gray average of data is m _to, gray variance is s _to, other moment t _fthe gray average of data is m _tf, gray variance is s _tf, gray-scale value is g _tf, then new time proximity sequence data corresponding grey scale value is:

$g_{\mathrm{new}}=\frac{[g_{\mathrm{tf}}-m_{\mathrm{tf}}]s_{\mathrm{to}}}{s_{\mathrm{tf}}}+m_{\mathrm{to}}.$

5. the time domain method for reconstructing of remote sensing sequence data according to claim 4, is characterized in that, the specific implementation of step 4 comprises following sub-step:

Step 4.1: to reconstruction moment t _odata are recovered;

First be the new time series g of step 3 gained _newin other time data distribute weights according to the related coefficient (CC) of step 2 gained, its weight computing mode is:

$w_i=\frac{{\mathrm{CC}}_i}{\underset{i=1}{\overset{q}{\mathrm{\Σ}}}{\mathrm{CC}}_i};$

Wherein, q is the number without missing data in time series, CC _iit is the related coefficient of i-th data and the public cloud-free area of data to be reconstructed;

Next converts data to the form that sine and cosine is added, and answers gray-scale value g to the new time proximity sequence pair of step 3 gained _newbe weighted and solve, its implementation is:

$\begin{array}{c}{\left|\right|g\left(t\right)-[a_0+\underset{j=1}{\overset{\mathrm{nf}}{\mathrm{\Σ}}}[a_j\cos \left(2\mathrm{\π}{f}_{j}t\right)+b_j\sin \left({2\mathrm{\πf}}_{j}t\right)\left]\right]\left|\right|}_w\\ =\∫{[g\left(t_i\right)-[a_0+\underset{j=1}{\overset{\mathrm{nf}}{\mathrm{\Σ}}}[a_j\cos \left(2\mathrm{\π}{f}_{j}t\right)+b_j\sin \left(2\mathrm{\π}{f}_{j}t\right)]\left]\right]}^2{w}_i\mathrm{dy}\end{array};$

Wherein, nf is harmonic wave number, a ₀for coefficient when harmonic wave number is 0, a _jand b _jfor the coefficient of trigonometric function each several part, f _jfor frequency, t _jfor timing node, y is time proximity sequence data;

Step 4.2: travel through all single picture point time sequences, repeating above-mentioned steps can rebuild remote sensing time series data, the reconstruction mark problem caused with the grey level's difference eliminated between different images.