CN102073037A - Iterative current inversion method based on adaptive threshold selection technique - Google Patents

Abstract

The present invention provides an iterative current inversion method based on an adaptive threshold selection technique, which comprises the following steps: (1) collecting radar images and carrying out A/D (analogue/digital) conversion to obtain 32 original radar images; (2) obtaining images within Cartesian coordinates and then obtaining an image spectrum through Fourier transformation; (3), determining initial estimates, only taking influence of 0-order secondary waves into consideration, and obtaining a rough flow through inversion; and (4) iterating fitting flow information by adopting a average weighted least squares method based on the adaptive threshold selection technique, taking influence of the 0-order secondary waves and 1-order secondary waves into consideration. By adopting the method, wave signals can be automatically identified, wave energy can be accurately separated from the image spectrum, the insufficient empirical value can be made up, the current trend can be better reflected, and the purposes of fewer errors and better stability can be achieved.

Description

Iteration ocean current inversion method based on the adaptive threshold selecting technology

Technical field

The present invention relates to a kind of remote sensing technology, particularly based on the ocean remote sensing technology of X-band navar.

Background technology

Flow measurement instrument (current meter, acoustics ocean current profile instrument) can be measured the stream information on surface, sea, but it is subjected to the destruction of artificial or disaster easily, so only use at present in the Yu Haiyang scientific investigation, China does not also form real-time and effective businessization monitoring capability to the flow field, sea at present.

X-band radar not only can monitor mobile target, and it is operated in short pulse pattern following time, and it can measure ocean current and wave information, and the measurement of flow velocity (size, direction) is the gordian technique of X-band radar wave monitoring.Compare with other sea sensor (wave buoy, flow measurement instrument), the X-band radar investigative range is wider, and is safe, safeguards advantages such as simple, can be used as the substitute of wave buoy and flow measurement instrument.Domestic existing how tame unit carries out Primary Study with regard to X-band radar flow measurement problem.At present, Chang Yong ocean current inversion technique is achieved as follows:

(1) (x, y t) arrive wave number frequency field F through 3 dimension discrete Fourier transformations for the time of radar image and spatial sequence g ⁽³⁾(k _x, k _y, ω), wherein, x, y represent two components in radar image space respectively, t represents the time shaft of radar image, k _x, k _yBe respectively the x and the y component of wave number, ω is the frequency of wave.

(2) because the wave number and the frequency of wave itself meet dispersion relation, as follows:

$S\left(k\right)=\sqrt{g\left|k\right|\tanh \left(\right|k\left|h\right)}+k\·u---\left(1\right)$

Wherein g is an acceleration of gravity, k=(k _x, k _y) be wave number, h is the depth of water, u is a flow velocity.

In the image spectrum, select an energy threshold C, require energy value to be generally dozens of greater than the number of the discrete point of C, use these conventional point,, use the least square fitting of energy weighting to obtain rough flow velocity according to the characteristic of wave number of wave own and frequency.

(3) at first according to the flow velocity that obtains, determine 0,1 rank dispersion relation.Rule of thumb choose energy threshold C then ₁, require C ₁Much smaller than C, choose greater than C ₁The spread in energy point, judge respectively that then it still is 1 order ripple that these points belong to 0 rank, simulate flow velocity according to the energy weighted least-squares method, new flow velocity is determined 0,1 new rank dispersion relation, repeat this step, will obtain constantly flow velocity accurately.

$w_p=(p+1)\sqrt{\frac{\mathrm{gk}}{p+1}\tanh \left(\frac{\mathrm{kh}}{p+1}\right)}+k\·u---\left(2\right)$

Wherein, p=0,1 o'clock respectively corresponding 0,1 rank dispersion relation.

This method mainly has the following disadvantages: 1. the extra large subsurface runoff that measures of X-band radar should be an average flow information in the zone of analyzing, and the stream that the least square method match of energy weighting obtains, and tends to the stream of the big ripple of energy; 2. rely on the method applicability of experience selected threshold bad, the empirical value difference of different waters, sea situation needs to obtain through a large amount of experiments, and versatility is bad.

Summary of the invention

It is a kind of accurate to fluid-velocity survey to the object of the present invention is to provide, and versatility is good, can improve the iteration ocean current inversion method based on the adaptive threshold selecting technology that X-band radar is measured the precision of extra large subsurface runoff.

The object of the present invention is achieved like this:

(1) collection of radar image and A D be converted to 32 original radar images,

The X-band navar is operated under the short pulse, the electromagnetic wave of emission and the capillary wave generation Bragg diffraction on sea, back scattering are received machine and receive, and vision signal is converted to digital signal through A/D, vision signal between two stem signals is formed a width of cloth radar image, select the continuous radar image g of 32 width of cloth (x, y, t), wherein, x, y are respectively two components in radar image space, and t is the time component;

(2) obtain the image under the Cartesian coordinates and it is obtained image spectrum as Fourier transform;

In radar image, select to need analyze rectangle frame, use the closest approach interpolation realize in the rectangle frame radar image by polar coordinates to the Cartesian coordinates conversion table be shown g1 (x, y, t); (x, y t) become the wave number frequency field to radar image sequence g1 in the rectangle frame after discrete 3 dimension Fourier transforms;

$\mathrm{\γ}(k,w)={\∫}_0^{L_X}{\∫}_0^{L_Y}{\∫}_0^{T}g(x,y,t)\exp \left[i(k_{x}x+k_{y}y-\mathrm{wt})\right]\mathrm{dxdydt}$

K=(k wherein _x, k _y), L _x, L _y, T is respectively the 180 degree fuzzy problems of considering actual conditions behind the length and width of rectangular area and the seasonal effect in time series total length of time three-dimensional Fourier transform and eliminating spectrum, only keeps the part of w＞0 so the energy spectrum of 3 d image and is:

$F^{\left(3\right)}(k,w)=\frac{1}{L_X{L}_{Y}T}{\mathrm{\γ}}^2(k,w);$

(3) influence of 0 order ripple is only considered in initial valuation, and inverting obtains rough stream;

In the wave number frequency field, choose the bigger point of dozens of energy value, according to the dispersion relation between Ocean-wave Signal wave number and frequency, use the average weighted least square method, simulate stream information, basic mathematical model is as follows:

In order to draw subsurface runoff u, with formula

${\mathrm{<figure-callout\; id="2"\; label="Q"\; filenames="HDA0000042754790000013.png,HDA0000042754790000021.png"\; state="\{\{state\}\}">Q</figure-callout>}}^2=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{({\mathrm{\&omega;}}_i-S\left(k_i\right))}^2$

To u _x, u _yAsk the single order partial derivative, and be zero:

$\&PartialD;{\mathrm{<figure-callout\; id="2"\; label="Q"\; filenames="HDA0000042754790000013.png,HDA0000042754790000021.png"\; state="\{\{state\}\}">Q</figure-callout>}}^2/\&PartialD;u_x=0$

$\&PartialD;{\mathrm{<figure-callout\; id="2"\; label="Q"\; filenames="HDA0000042754790000013.png,HDA0000042754790000021.png"\; state="\{\{state\}\}">Q</figure-callout>}}^2/\&PartialD;u_y=0$

Following formula is write as the form of matrix:

$\left(\begin{array}{cc}{D}_{\mathrm{xx}}& D_{\mathrm{xy}}\\ D_{\mathrm{yx}}& D_{\mathrm{yy}}\end{array}\right)\left(\begin{array}{c}{u}_x\\ u_y\end{array}\right)=\left(\begin{array}{c}{b}_x\\ b_y\end{array}\right)$

Be abbreviated as:

Du＝b

The element of matrix D is

$D_{\mathrm{xx}}=\underset{i=1}{\overset{\mathrm{<figure-callout\; id="2"\; label="n\; k\; ix"\; filenames="HDA0000042754790000013.png,HDA0000042754790000021.png"\; state="\{\{state\}\}">n</figure-callout>}}{\mathrm{\&Sigma;}}}{k}_{\mathrm{ix}}^{}{}_2^{},D_{\mathrm{xy}}=D_{\mathrm{yx}}=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{k}_{\mathrm{ix}}{k}_{\mathrm{iy}}$

$D_{\mathrm{yy}}=\underset{i=1}{\overset{\mathrm{<figure-callout\; id="2"\; label="n\; k\; iy"\; filenames="HDA0000042754790000013.png,HDA0000042754790000021.png"\; state="\{\{state\}\}">n</figure-callout>}}{\mathrm{\&Sigma;}}}{k}_{\mathrm{iy}}^{}{}_2^{}$

The element of vector b is

$b_x=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{k}_{\mathrm{ix}}{w}_{\mathrm{iD}},$ $b_y=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{k}_{\mathrm{iy}}{w}_{\mathrm{iD}}$

Wherein, w _ID=w _i-ζ (k _i) in order to calculate stream u, it is D that D is necessary for nonsingular matrix (det (D) ≠ 0) ^-1Exist, at this moment

u＝D ^-1·b；

(4) based on the adaptive threshold selecting technology, consider the influence of 0,1 order ripple, use average weighted least square method iterative fitting stream information;

1. use the adaptive threshold selecting technology, obtain threshold value C _It, greater than C _ItEnergy just comprised 0 order and 1 order Wave energy, number is N ₁, the stream that utilizes initial valuation to obtain is by formula

Calculate 0 order and 1 order wave frequency;

2. judge this N that surveys ₁Individual energy point meets 0 order or 1 order wave dispersion relation, if | w _i-w (k _i) |＜| w _i-w ₁(k _i) |, then this frequency meets 0 order ripple; If | w _i-w (k _i) |＞| w _i-w ₁(k _i) |, then this frequency meets 1 order ripple, with judge good data based different minimal value function suc as formula Use average least square method, match obtains new subsurface runoff;

3. use new subsurface runoff and make up bandpass filter, obtain new wave signal to noise ratio (S/N ratio), use the adaptive threshold selecting technology, obtain new threshold value;

4. the subsurface runoff that will newly obtain, the substitution formula

Obtain 0 new order and 1 order wave frequency, repeat above-mentioned steps, will obtain constantly accurate subsurface runoff.

Described adaptive threshold selecting technology is:

1. in the short period of time, the variation of wave on time and space is less, and the signal to noise ratio (S/N ratio) of a preceding sequence chart is the signal to noise ratio snr of this sequence chart in approximate thinking;

2. the image that is obtained by 3 dimension Fourier transforms is composed F ⁽³⁾(k, w) and signal to noise ratio snr obtain the Wave energy P that contains in the image spectrum

$P=\frac{\mathrm{SNR}+1}{\mathrm{SNR}}\underset{i}{\overset{N_{\mathrm{kx}}}{\mathrm{\&Sigma;}}}\underset{j}{\overset{N_{\mathrm{ky}}}{\mathrm{\&Sigma;}}}\underset{k}{\overset{N_w}{\mathrm{\&Sigma;}}}{F}^{\left(3\right)}(k_{\mathrm{ix}},k_{\mathrm{jy}},w_k)\&CenterDot;{\mathrm{\&Delta;k}}_x{\mathrm{\&Delta;k}}_y\mathrm{\&Delta;w}$

3. energy points all in the image spectrum is sorted according to from big to small order, the energy of ordering back i sampled point is P _i,, therefore can make because Wave energy point is generally big than the ground unrest energy point

$P=\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}{P}_i$

Obtain the energy point number m that is used for LSM by following formula, the energy value of m sampled point is P _m, P then _mBe iteration energy threshold C _It

Main technical points of the present invention is:

(1) initial valuation obtains rough flow velocity (size, direction).Suppose that the wave space is even in the rectangle frame, the time is stable, and then the wave model meets Gaussian distribution.When not having under the flow rate conditions, the dispersion relation of extra large capillary wave and gravity wave is

$\mathrm{\&zeta;}\left(k\right)\&ap;\sqrt{g\left|k\right|\tanh \left(\right|k\left|h\right)}---\left(3\right)$

Wherein, ζ (k) is a wave frequency, and k is a wave-number vector, and h is the depth of water, and g is an acceleration of gravity.If a subsurface runoff u with respect to radar is arranged, just introduced the Doppler shift item in the frequency, dispersion relation becomes:

S(k)＝ζ+k·u＝ζ+|k||u|cosθ(4)

Wherein, S (k) is theoretical wave frequency, and ku is a Doppler shift, u=(u _x, u _y) be the relative velocity between wave field and radar antenna platform, comprise Texas tower movement velocity (such as ship's speed) and subsurface runoff speed vector, θ is an angle between wave number k and the u, when θ=90, the Doppler shift item is 0, have only stream u could influence Doppler shift when important on the direction of wave number k, when Texas tower was static, u referred to extra large subsurface runoff.

Can obtain rough stream information by the average weighted least square method.The minimal value function definition of average weighted least square method is as the formula (5):

${\mathrm{<figure-callout\; id="2"\; label="Q"\; filenames="HDA0000042754790000013.png,HDA0000042754790000021.png"\; state="\{\{state\}\}">Q</figure-callout>}}^2=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{({\mathrm{\&omega;}}_i-S\left(k_i\right))}^2---\left(5\right)$

Wherein, ω _iBe the wave frequency of radar observation, n is the energy point number that is used for least square method.

(2) the automatic selecting technology of threshold value obtains threshold value C _It

Big when wind speed, when wave height is higher, the echo of radar image is stronger, the signal to noise ratio (S/N ratio) of Ocean-wave Signal is bigger, the quality of the data of this moment is better.Otherwise the Ocean-wave Signal signal to noise ratio (S/N ratio) is less, and the quality of data is poor.The quality of the signal to noise ratio (S/N ratio) judgment data quality by Ocean-wave Signal, selected threshold C according to this _It

At first provide the definition of wave image energy signal to noise ratio (S/N ratio):

$\mathrm{SNR}=\frac{\mathrm{SIG}}{\mathrm{BGN}}---\left(6\right)$

$\mathrm{SIG}=\underset{i}{\overset{N_{\mathrm{kx}}}{\mathrm{\&Sigma;}}}\underset{j}{\overset{N_{\mathrm{ky}}}{\mathrm{\&Sigma;}}}{F}^{\left(2\right)}(k_{\mathrm{ix}},k_{\mathrm{jy}}){\mathrm{\&Delta;k}}_x{\mathrm{\&Delta;k}}_y---\left(7\right)$

$\mathrm{BGN}=\underset{i}{\overset{N_{\mathrm{kx}}}{\mathrm{\&Sigma;}}}\underset{j}{\overset{N_{\mathrm{ky}}}{\mathrm{\&Sigma;}}}\underset{k}{\overset{N_{\mathrm{\&omega;}}}{\mathrm{\&Sigma;}}}{F}^{\left(3\right)}(k_{\mathrm{ix}},k_{\mathrm{jy}},{\mathrm{\&omega;}}_k){\mathrm{\&Delta;k}}_x{\mathrm{\&Delta;k}}_y\mathrm{\&Delta;\&omega;}---\left(8\right)$

$-\underset{i}{\overset{N_{\mathrm{kx}}}{\mathrm{\&Sigma;}}}\underset{j}{\overset{N_{\mathrm{ky}}}{\mathrm{\&Sigma;}}}{F}^{\left(2\right)}(k_{\mathrm{ix}},k_{\mathrm{jy}}){\mathrm{\&Delta;k}}_x{\mathrm{\&Delta;k}}_y$

$F^{\left(2\right)}\left(k\right)=2{\&Integral;}_{\mathrm{\&omega;}>0}{F}^{\left(3\right)}(k_x,k_y,\mathrm{\&omega;})\mathrm{\&delta;}(\mathrm{\&omega;}-{\mathrm{\&omega;}}_0)\mathrm{d\&omega;}---\left(9\right)$

Wherein, SIG is the wave spectrum energy, and BGN is a ground unrest, F ⁽³⁾Be the 3 d image spectrum, F ⁽²⁾(k) be 2 dimension image spectrums, N _Kx, N _Ky, N _ωBe the scope of spectrum, Δ k _x, Δ k _yBe wavenumber resolution, Δ ω is a frequency resolution, δ (ω-ω ₀) be bandpass filter.

The basic realization of threshold value automatically selecting method, 7 explanations in conjunction with the accompanying drawings.

1. in the short period of time, it promptly is stable in time that the signal to noise ratio (S/N ratio) of wave changes less, and therefore the signal to noise ratio (S/N ratio) that obtains before the hypothesis is the signal to noise ratio snr of this sequence chart.

2. the image that is obtained by 3 dimension Fourier transforms is composed F ⁽³⁾(k, w) and signal to noise ratio snr can obtain the Wave energy P that contains in the image spectrum.

$P=\frac{\mathrm{SNR}+1}{\mathrm{SNR}}\underset{i}{\overset{N_{\mathrm{kx}}}{\mathrm{\&Sigma;}}}\underset{j}{\overset{N_{\mathrm{ky}}}{\mathrm{\&Sigma;}}}\underset{k}{\overset{N_w}{\mathrm{\&Sigma;}}}{F}^{\left(3\right)}(k_{\mathrm{ix}},k_{\mathrm{jy}},w_k)\&CenterDot;{\mathrm{\&Delta;k}}_x{\mathrm{\&Delta;k}}_y\mathrm{\&Delta;w}---\left(10\right)$

3. energy points all in the image spectrum is sorted according to from big to small order, the energy of ordering back i sampled point is P _i,, therefore can make because Wave energy point is generally big than the ground unrest energy point

$P=\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}{P}_i---\left(11\right)$

Can obtain the energy point number m that is used for LSM by following formula, the energy value of m sampled point is P _m, P then _mBe iteration energy threshold C _It

(3) based on the iteration valuation of the automatic selecting technology of threshold value.

At this, we only consider the situation of 0,1 order harmonic wave.

1. application self-adapting selection of threshold technology, selected threshold C _iT can think greater than C _ItEnergy partly comprised 0 order and 1 order Wave energy, number is N ₁The stream that utilizes initial valuation to obtain calculates 0 order and 1 order wave frequency by formula (12) respectively.

$w_p=(p+1)\sqrt{\frac{\mathrm{gk}}{p+1}\tanh \left(\frac{\mathrm{kh}}{p+1}\right)}+k\&CenterDot;u---\left(12\right)$

The dispersion relation of p=0,1 respectively corresponding 0,1 order ripple.

2. judge this N that actual measurement obtains ₁Individual energy point belongs to 0 order or 1 order ripple, if | w _i-w (k _i) |＜| w _i-w ₁(k _i) |, then this frequency meets 0 order ripple; If | w _i-w (k _i) |＞| w _i-w ₁(k _i) |, then this frequency meets 1 order ripple.With judging good data based different minimal value function, use least square method, can obtain new subsurface runoff.

3. use new subsurface runoff and make up bandpass filter, obtain new Ocean-wave Signal signal to noise ratio (S/N ratio), application self-adapting selection of threshold technology obtains new threshold value.

4. the subsurface runoff that will newly obtain in the substitution formula (12), obtains 0 new order and 1 order wave dispersion relation, repeats above-mentioned steps, will obtain constantly accurate subsurface runoff.

The stream inverting new method advantage of choosing based on automatic threshold compared with prior art of invention:

(1) different wind speed and zone down in the radar sequence chart Ocean-wave Signal be different, mainly be to embody by signal to noise ratio (S/N ratio).The energy of wave in the image sequence can be calculated by the gross energy of signal to noise ratio (S/N ratio) and image, threshold value can be drawn thus.The present invention proposes the automatic selecting technology of threshold value according to above thinking, it can discern the Ocean-wave Signal quality automatically, and by iterative technique, isolates Wave energy accurately from the image spectrum, has remedied the deficiency of empirical value.

(2) as can be known, from Fig. 1,2, table 1 when flow velocity is big (＞0.3m/s), the more former methods and results of stream inversion result of the inventive method has two aspect advantages: the trend that 1. can better reflect stream; 2. has better stability of littler sum of errors.

(3) flow velocity hour (＜0.2m/s), signal a little less than, extracting stream information from weak signal is a difficult problem.The more former method of the inversion result of method of the present invention has been significantly improved shown in Fig. 3,4, table 1.

Description of drawings

When Fig. 1 flow velocity was big, method, former method and spot sensor flow rate result were relatively among the present invention.

When Fig. 2 flow velocity was big, method, former method and spot sensor flowed to the result relatively among the present invention.

Fig. 3 flow velocity hour, method, former method and the comparison of spot sensor flow rate result among the present invention.

Fig. 4 flow velocity hour, method, former method and spot sensor flow to result's comparison among the present invention.

Fig. 5 technology implementation scheme of the present invention process flow diagram.

Fig. 6 is based on the iteration stream inverting process flow diagram of adaptive threshold selecting technology.

Fig. 7 adaptive threshold selecting technology implementing procedure figure.

When Fig. 8 table 1 flow velocity is big, each methods and results statistical study.

Fig. 9 table 2 flow velocity hour, each methods and results statistical study.

Embodiment

Below in conjunction with Fig. 5,6,7, the present invention is done more detailed introduction:

(5) collection of radar image and A D be converted to 32 original radar images.

The X-band navar is operated under the short pulse, and the capillary wave generation Bragg diffraction on the electromagnetic wave of emission and sea, back scattering are received machine and receive, vision signal is converted to digital signal through A/D, vision signal between two stem signals is formed a width of cloth radar image, in the present invention, needs to use the continuous radar image g (x of 32 width of cloth, y, t), wherein, x, y is respectively two components in radar image space, and t is the time component.

(6) obtain the image under the Cartesian coordinates and it is obtained image spectrum as Fourier transform.

In radar image, select to need analyze rectangle frame, use the closest approach interpolation realize in the rectangle frame radar image by polar coordinates to the Cartesian coordinates conversion table be shown g1 (x, y, t); (x, y t) become the wave number frequency field to radar image sequence g1 in the rectangle frame after discrete 3 dimension Fourier transforms.

$\mathrm{\&gamma;}(k,w)={\&Integral;}_0^{L_X}{\&Integral;}_0^{L_Y}{\&Integral;}_0^{T}g(x,y,t)\exp \left[i(k_{x}x+k_{y}y-\mathrm{wt})\right]\mathrm{dxdydt}---\left(13\right)$

K=(k wherein _x, k _y), L _x, L _y, T is respectively the 180 degree fuzzy problems of considering actual conditions behind the length and width of rectangular area and the seasonal effect in time series total length of time three-dimensional Fourier transform and eliminating spectrum, and we have only kept the part of w＞0 so the energy spectrum of 3 d image is:

$F^{\left(3\right)}(k,w)=\frac{1}{L_X{L}_{Y}T}{\mathrm{\&gamma;}}^2(k,w)---\left(14\right)$

(7) influence of 0 order ripple is only considered in initial valuation, and inverting obtains rough stream.

In the wave number frequency field, choose the bigger point of dozens of energy value, according to the dispersion relation between Ocean-wave Signal wave number and frequency, use the average weighted least square method, simulate stream information, basic mathematical model is as follows:

In order to draw subsurface runoff u, with formula (15)

${\mathrm{<figure-callout\; id="2"\; label="Q"\; filenames="HDA0000042754790000013.png,HDA0000042754790000021.png"\; state="\{\{state\}\}">Q</figure-callout>}}^2=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{({\mathrm{\&omega;}}_i-S\left(k_i\right))}^2---\left(15\right)$

To u _x, u _yAsk the single order partial derivative, and be zero:

$\begin{array}{c}\&PartialD;{\mathrm{<figure-callout\; id="2"\; label="Q"\; filenames="HDA0000042754790000013.png,HDA0000042754790000021.png"\; state="\{\{state\}\}">Q</figure-callout>}}^2/\&PartialD;u_x=0\\ \&PartialD;{\mathrm{<figure-callout\; id="2"\; label="Q"\; filenames="HDA0000042754790000013.png,HDA0000042754790000021.png"\; state="\{\{state\}\}">Q</figure-callout>}}^2/\&PartialD;u_y=0\end{array}---\left(16\right)$

Following formula is write as the form of matrix:

$\left(\begin{array}{cc}{D}_{\mathrm{xx}}& D_{\mathrm{xy}}\\ D_{\mathrm{yx}}& D_{\mathrm{yy}}\end{array}\right)\left(\begin{array}{c}{u}_x\\ u_y\end{array}\right)=\left(\begin{array}{c}{b}_x\\ b_y\end{array}\right)---\left(17\right)$

Be abbreviated as:

Du＝b (18)

The element of matrix D is

$D_{\mathrm{xx}}=\underset{i=1}{\overset{\mathrm{<figure-callout\; id="2"\; label="n\; k\; ix"\; filenames="HDA0000042754790000013.png,HDA0000042754790000021.png"\; state="\{\{state\}\}">n</figure-callout>}}{\mathrm{\&Sigma;}}}{k}_{\mathrm{ix}}^{}{}_2^{},D_{\mathrm{xy}}=D_{\mathrm{yx}}=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{k}_{\mathrm{ix}}{k}_{\mathrm{iy}}---\left(19\right)$

$D_{\mathrm{yy}}=\underset{i=1}{\overset{\mathrm{<figure-callout\; id="2"\; label="n\; k\; iy"\; filenames="HDA0000042754790000013.png,HDA0000042754790000021.png"\; state="\{\{state\}\}">n</figure-callout>}}{\mathrm{\&Sigma;}}}{k}_{\mathrm{iy}}^{}{}_2^{}$

The element of vector b is

$b_x=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{k}_{\mathrm{ix}}{w}_{\mathrm{iD}},b_y=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{k}_{\mathrm{iy}}{w}_{\mathrm{iD}}---\left(20\right)$

Wherein, w _ID=w _i-ζ (k _i) in order to calculate stream u, it is D that D is necessary for nonsingular matrix (det (D) ≠ 0) ^-1Exist, at this moment

u＝D ^-1·b (21)

(8) based on the adaptive threshold selecting technology, consider the influence of 0,1 order ripple, use average weighted least square method iterative fitting stream information.6 describe in conjunction with the accompanying drawings.The concrete enforcement of adaptive threshold selecting technology sees (5) part for details.

The flow velocity (size, direction) that obtains by the initial valuation of (3) step is inaccurate, needs to consider the influence of single order subwave, uses the inversion accuracy of the method raising stream of iteration, and its implementation 6 is done following explanation in conjunction with the accompanying drawings:

1. use the adaptive threshold selecting technology, obtain threshold value C _It, greater than C _ItEnergy just comprised 0 order and 1 order Wave energy, number is N ₁The stream that utilizes initial valuation to obtain calculates 0 order and 1 order wave frequency by formula (12).

2. judge this N that surveys ₁Individual energy point meets 0 order or 1 order wave dispersion relation, if | w _i-w (K _i) |＜| w _i-w ₁(k _i) |, then this frequency meets 0 order ripple; If | w _i-w (k _i) |＞| w _i-w ₁(k _i) |, then this frequency meets 1 order ripple.The data based different minimal value function that judgement is good is used average least square method suc as formula (15), can match obtain new subsurface runoff.

3. use new subsurface runoff and make up bandpass filter, obtain new wave signal to noise ratio (S/N ratio), use the adaptive threshold selecting technology, obtain new threshold value.

4. the subsurface runoff that will newly obtain, substitution formula (12) obtains 0 new order and 1 order wave frequency, repeats above-mentioned steps, will obtain constantly accurate subsurface runoff.

(9) in last step, use the adaptive threshold selecting technology, now 7 describe in conjunction with the accompanying drawings,

1. in the short period of time, the variation of wave on time and space is less, and what therefore can be similar to thinks that the signal to noise ratio (S/N ratio) of a preceding sequence chart is the signal to noise ratio snr of this sequence chart;

2. the image that is obtained by 3 dimension Fourier transforms is composed F ⁽³⁾(k, w) and signal to noise ratio snr can obtain the Wave energy P that contains in the image spectrum

$P=\frac{\mathrm{SNR}+1}{\mathrm{SNR}}\underset{i}{\overset{N_{\mathrm{kx}}}{\mathrm{\&Sigma;}}}\underset{j}{\overset{N_{\mathrm{ky}}}{\mathrm{\&Sigma;}}}\underset{k}{\overset{N_w}{\mathrm{\&Sigma;}}}{F}^{\left(3\right)}(k_{\mathrm{ix}},k_{\mathrm{jy}},w_k)\&CenterDot;{\mathrm{\&Delta;k}}_x{\mathrm{\&Delta;k}}_y\mathrm{\&Delta;w}---\left(22\right)$

3. energy points all in the image spectrum is sorted according to from big to small order, the energy of ordering back i sampled point is P _i,, therefore can make because Wave energy point is generally big than the ground unrest energy point

$P=\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}{P}_i---\left(23\right)$

Can obtain the energy point number m that is used for LSM by following formula, the energy value of m sampled point is P _m, P then _mBe iteration energy threshold C _It

The stream inversion method that the present invention is proposed is applied in the measured data of flow velocity when big (time period be the morning on the 06th October in 2009 10:36:00 to 16:40:00 in afternoon), and the result as shown in Figure 1, 2, statistic analysis result is as shown in table 1, and the more former method of the inventive method result has higher precision and stability as can be seen.Fig. 3,4 has provided hour result of the inventive method inverting stream (time period be the morning on the 06th October in 2009 08:16:00 to 09:20:00), statistics analysis in table 2, and the more former methods and results precision of the inventive method result improves a lot as can be known.

Claims (3)

1. iteration ocean current inversion method based on the adaptive threshold selecting technology is characterized in that:

(1) collection of radar image and A D be converted to 32 original radar images;

(2) obtain the image under the Cartesian coordinates and it is obtained image spectrum as Fourier transform;

(3) influence of 0 order ripple is only considered in initial valuation, and inverting obtains rough stream;

(4) based on the adaptive threshold selecting technology, consider the influence of 0,1 order ripple, use average weighted least square method iterative fitting stream information.

2. the iteration ocean current inversion method based on the adaptive threshold selecting technology according to claim 1 is characterized in that the concrete grammar of each step is:

(1) collection of radar image and A D be converted to 32 original radar images;

The X-band navar is operated under the short pulse, the electromagnetic wave of emission and the capillary wave generation Bragg diffraction on sea, back scattering are received machine and receive, and vision signal is converted to digital signal through A/D, vision signal between two stem signals is formed a width of cloth radar image, select the continuous radar image g of 32 width of cloth (x, y, t), wherein, x, y are respectively two components in radar image space, and t is the time component;

(2) obtain the image under the Cartesian coordinates and it is obtained image spectrum as Fourier transform;

In radar image, select to need analyze rectangle frame, use the closest approach interpolation realize in the rectangle frame radar image by polar coordinates to the Cartesian coordinates conversion table be shown g1 (x, y, t); (x, y t) become the wave number frequency field to radar image sequence g1 in the rectangle frame after discrete 3 dimension Fourier transforms;

$\mathrm{\&gamma;}(k,w)={\&Integral;}_0^{L_X}{\&Integral;}_0^{L_Y}{\&Integral;}_0^{T}g(x,y,t)\exp \left[i(k_{x}x+k_{y}y-\mathrm{wt})\right]\mathrm{dxdydt}$

K=(k wherein _x, k _y), L _x, L _y, T is respectively the 180 degree fuzzy problems of considering actual conditions behind the length and width of rectangular area and the seasonal effect in time series total length of time three-dimensional Fourier transform and eliminating spectrum, only keeps the part of w＞0 so the energy spectrum of 3 d image and is:

$F^{\left(3\right)}(k,w)=\frac{1}{L_X{L}_{Y}T}{\mathrm{\&gamma;}}^2(k,w);$

(3) influence of 0 order ripple is only considered in initial valuation, and inverting obtains rough stream;

In the wave number frequency field, choose the bigger point of dozens of energy value, according to the dispersion relation between Ocean-wave Signal wave number and frequency, use the average weighted least square method, simulate stream information, basic mathematical model is as follows:

In order to draw subsurface runoff u, with formula

$Q^2=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{({\mathrm{\&omega;}}_i-S\left(k_i\right))}^2$

To u _x, u _yAsk the single order partial derivative, and be zero:

$\begin{array}{c}\&PartialD;Q^2/\&PartialD;u_x=0\\ \&PartialD;Q^2/\&PartialD;u_y=0\end{array}$

Following formula is write as the form of matrix:

$\left(\begin{array}{cc}{D}_{\mathrm{xx}}& D_{\mathrm{xy}}\\ D_{\mathrm{yx}}& D_{\mathrm{yy}}\end{array}\right)\left(\begin{array}{c}{u}_x\\ u_y\end{array}\right)=\left(\begin{array}{c}{b}_x\\ b_y\end{array}\right)$

Be abbreviated as:

Du＝b

The element of matrix D is

$D_{\mathrm{xx}}=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{k}_{\mathrm{ix}}^2,D_{\mathrm{xy}}=D_{\mathrm{yx}}=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{k}_{\mathrm{ix}}{k}_{\mathrm{iy}}$

$D_{\mathrm{yy}}=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{k}_{\mathrm{iy}}^2$

The element of vector b is

$b_x=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{k}_{\mathrm{ix}}{w}_{\mathrm{iD}},b_y=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{k}_{\mathrm{iy}}{w}_{\mathrm{iD}}$

Wherein, w _ID=w _i-ζ (k _i) in order to calculate stream u, it is D that D is necessary for nonsingular matrix (det (D) ≠ 0) ^-1Exist, at this moment

u＝D ^-1·b；

(4) based on the adaptive threshold selecting technology, consider the influence of 0,1 order ripple, use average weighted least square method iterative fitting stream information;

1. use the adaptive threshold selecting technology, obtain threshold value C _It, greater than C _ItEnergy just comprised 0 order and 1 order Wave energy, number is N ₁, the stream that utilizes initial valuation to obtain is by formula

Calculate 0 order and 1 order wave frequency;

2. judge this N that surveys ₁Individual energy point meets 0 order or 1 order wave dispersion relation, if | w _i-w (k _i) |＜| w _i-w ₁(k _i) |, then this frequency meets 0 order ripple; If | w _i-w (k _i) |＞| w _i-w ₁(k _i) |, then this frequency meets 1 order ripple, with judge good data based different minimal value function suc as formula

Use average least square method, match obtains new subsurface runoff;

3. use new subsurface runoff and make up bandpass filter, obtain new wave signal to noise ratio (S/N ratio), use the adaptive threshold selecting technology, obtain new threshold value;

4. the subsurface runoff that will newly obtain, the substitution formula

Obtain 0 new order and 1 order wave frequency, repeat above-mentioned steps, will obtain constantly accurate subsurface runoff.

3. the iteration ocean current inversion method based on the adaptive threshold selecting technology according to claim 2 is characterized in that described adaptive threshold selecting technology is:

1. in the short period of time, the variation of wave on time and space is less, and the signal to noise ratio (S/N ratio) of a preceding sequence chart is the signal to noise ratio snr of this sequence chart in approximate thinking;

2. the image that is obtained by 3 dimension Fourier transforms is composed F ⁽³⁾(k, w) and signal to noise ratio snr obtain the Wave energy P that contains in the image spectrum

$P=\frac{\mathrm{SNR}+1}{\mathrm{SNR}}\underset{i}{\overset{N_{\mathrm{kx}}}{\mathrm{\&Sigma;}}}\underset{j}{\overset{N_{\mathrm{ky}}}{\mathrm{\&Sigma;}}}\underset{k}{\overset{N_w}{\mathrm{\&Sigma;}}}{F}^{\left(3\right)}(k_{\mathrm{ix}},k_{\mathrm{jy}},w_k)\&CenterDot;{\mathrm{\&Delta;k}}_x{\mathrm{\&Delta;k}}_y\mathrm{\&Delta;w}$

3. energy points all in the image spectrum is sorted according to from big to small order, the energy of ordering back i sampled point is P _i,, therefore can make because Wave energy point is generally big than the ground unrest energy point

$P=\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}{P}_i$

Obtain the energy point number m that is used for LSM by following formula, the energy value of m sampled point is P _m, P then _mBe iteration energy threshold C _It

Images