CN100349183C - Two-dimensional image data separation with multiple scale information - Google Patents

Abstract

The present invention relates to a separation method of two-dimensional image data with the information of multiple scales. A two-dimensional image which comprises the information of various characteristic scales is divided into a plurality of 'layers' with separate scales, wherein each of the 'layers' comprises only the information of a goal object with a relatively single scale, so the tranquilized pretreatment of the two-dimensional data and the separation of superimposed information of the different scales are realized, and data source pretreatment is provided to the mutual comparison and the research among the phenomena of the scales in following data application phase. The core technical scheme of the separation method of the data comprises two points: 1. the separated image is composed of two-dimensional data, and both the image which needs to be separated and the separated images are composed of the two-dimensional data; 2. the images comprise the information of multiple scales. The present invention has the aim that one two-dimensional image is divided into a plurality of 'layers' with separate scales, and each of the 'layers' comprises only the information of a goal object with a relatively single scale.

Description

Two-dimensional image data separation method with multi-scale information

Technical field

The present invention relates to a kind ofly carry out the method for data separating, specifically the view data of two dimension is separated by each scale calibration at the formed two dimensional image of electromagnetism echo.

Background technology

In fields such as ocean expedition and geologic prospecting, adopt electromagnetism echo remote sensing technology at present comparatively at large.The electromagnetism echo-wave imaging is the electromagnetic radiation characteristic that utilizes the object uniqueness, when photoelectricity or microwave detector receive the electromagnetic wave of object emission, reflection or stimulated emission, can convert electromagnetic radiation to electric signal, electric signal can be exchanged into image through further handling.The electromagnetic wave that each part (being called picture dot) of detection thing is produced has directly reflected the physical characteristics of surveying thing, the signal that is converted is reflected on the two dimensional image, features such as the intensity of image picture elements value, distribution are directly corresponding with the character of surveying thing, can reflect the characteristic information of object.Because the picture dot character of a certain scene or object has the otherness that is different from other, thereby causes its pixel value that produces also to have otherness, the randomness of these numerical value is very strong, so view data has not stationarity.

Because of object often has different characteristic dimensions according to the different criteria for classifying, such as when studying oceanographic phenomena, usually by its space scale and physical connotation it is decomposed into small scale (100m) ocean wave motion, than large scale (1km-10km) ocean current field and the inhomogeneous part of wind field, and three parts such as large scale mean field and circulation field.Utilize the stack of the oceanographic phenomena of the above-mentioned often three kinds of different scale information of image that the electromagnetism echo technique obtains.Be reflected on the image during information of same of different scale, be difficult to directly from image, obtain the marine physics information of different scale.When one of them yardstick phenomenon of research, will inevitably be subjected to the interference of other yardstick information.In like manner, in a lot of fields, the information of utilizing electromagnetism echo technique detection object to obtain generally all shows in same the two dimensional image, the target information of the same characteristic dimension that object has just is difficult to express singlely, also just must have influence on the data application stage of the follow-up flow process of data extract, thereby have influence on the positive stop of target detection thing and predict the result.

Existing a kind of data separating disposal route that is applied to above-mentioned field is the Empirical Mode Decomposition method (Empirical ModeDecomposition is hereinafter to be referred as EMD) that was proposed by Dr.Norden Huang in 1998.This method is the non-stationary signal sequence with a complexity, resolves into the stack of several steady data Layers and trend term step by step.Wherein each steady data Layer is called as a natural mode function (intrinsic modefunction), and a steady data Layer is represented the dimensional properties of a data-signal.Thereby can directly extract the image information based on various yardsticks of data-signal inside by this kind one-dimensional facture.

The shortcoming of said method is only to handle one dimension non-stationary information, promptly mainly is to handle at the actual measurement sequence of points data of object.Also only can just can handle after with its generate one-dimensional data sequence 2-D data, this has obviously ignored the space two-dimensional correlation properties of each information word.And in electromagnetic surveying applications such as ocean expedition, geologic prospecting, medical diagnosiss, the object of processing is often based on two-dimensional image data, and the 3 D stereo target also is processed into the 2-D data of some sections usually and studies.The target that institute's desire is separated, as use synthetic-aperture radar (SAR, Synthetic Aperture Radar) the ocean remote sensing image of Sheng Chenging is surveyed submarine topography information, no matter from form still with the correlativity of peripheral information, all have two-dimensional characteristics, therefore, only handle the electromagnetism echo image data of two-dimentional non-stationary with the one dimension separation method, obviously can't meet the demands the information that the accuracy of its result is relatively poor, tend to be doped with other characteristic dimensions.

Summary of the invention

Two-dimensional image data separation method with multi-scale information of the present invention, its purpose is to address the above problem and is not enough and design has the decomposition method that embodies two-dimensional image information space distribution state to handle 2-D data, thereby a two dimensional image that includes various features yardstick information can be separated into " the figure layer " that several yardsticks separate, wherein each opens the destination object information that " figure layer " only contains single relatively yardstick, thereby realize separating of 2-D data tranquilization pre-service and different scale overlapped information, for mutual comparison of each yardstick phenomenon and the research of follow-up data application stage provides the data source pre-service.

For echo image data the processing stage can be divided into that data are obtained, information extraction and data use this three phases, described two-dimensional image data separation method with multi-scale information is applied to the information extraction stage.At the two dimensional image that data acquisition phase feeds back to, therefrom isolate the view data of various yardsticks.

Described two-dimensional image data separation method with multi-scale information, its core technology scheme include following 2 points:

The one, the image that is separated is made of 2-D data, the separating treatment of and what separate all is images that 2-D data constituted.

The 2nd, include multi-scale information in the image, goal of the invention is separated into several yardsticks with a two dimensional image and separates " figure layer ", and wherein each opens the destination object information that " figure layer " only contains single relatively yardstick.

Based on above-mentioned 2 points, use the single figure layer that described two-dimensional image data separation method is obtained, owing to do not comprise the stack of other different scale information, so should be a data smooth change image.

When the separate picture that obtains approached single figure layer, the difference between two images was more little, and its similarity is big more.

Set a similarity coefficient PJ at this and express similarity between the above-mentioned image.By theoretical or a large amount of experiments predicts data, can obtain the yardstick information characteristic of generally a certain scene or object, so similarity coefficient PJ can set an area boundary in advance, the area boundary parameter setting is α ₁And α ₂

Two-dimensional image data separation method with multi-scale information of the present invention is according to the area boundary parameter alpha of setting ₁And α ₂, the single figure layer that similarity coefficient PJ falls in the area boundary scope is separated.

Described two-dimensional image data separation method with multi-scale information, implementation step is:

Obtain two-dimensional image data, promptly obtain in the image original two dimensional data matrix F (x, y).

Determine above-mentioned two-dimensional image data matrix F (x, y) each Local Extremum, and connect by certain syntople each subdivision triangle is carried out interpolation to guarantee second-order smooth on the border, carry out surface fitting for maximum point in the image and minimum point, thereby form coenvelope face Xmax (x, y) and lower envelope face Xmin (x y), obtains an average curved surface m by upper and lower envelope surface ₁(x, y);

(x gets rid of average curved surface m in y) from original image F _k(x y) obtains error image h _1k, work as k=1, promptly by average curved surface m ₁(x y) obtains error image h ₁₁

With error image h ₁₁Be the original image of data separating next time, repeat above-mentioned detachment process and obtain average curved surface m ₂(x is y) with error image h ₁₂, obtain average curved surface m until repeating above-mentioned separation method k time _k(x is y) with error image h _1k

Wherein, the k secondary data is separated the error image h that obtains _1kBe an above error image h _{1 (k-1)}For original image, wherein k is constant and k 〉=2.

As difference image h _1kThe average curved surface m that obtains with the last time _K-1(x, y) determined error image h _{1 (k-1)}, when satisfying following similarity coefficient PJ function expression, above-mentioned data separating process end:

$\mathrm{PJ}=\underset{x=1}{\overset{W}{\mathrm{\Σ}}}\underset{y=1}{\overset{H}{\mathrm{\Σ}}}\left[\frac{{|(h_{m(k-1)}(x,y)-h_{\mathrm{mk}}(x,y)|}^2}{h_{m(k-1)}^2(x,y)}\right]\∈[{\mathrm{\α}}_1,{\mathrm{\α}}_2]---\left(3.3\right)$

Wherein, m is the number of plies when pre-treatment;

K is the number of times of this layer re-treatment, and k 〉=2;

W and H are respectively picture traverse and height;

And α ₁And α ₂Be the steady data Layer area boundary parameter of single yardstick information given before the data separating, α ₁And α ₂Be constant, and α ₁＜α ₂

By above-mentioned two-dimensional image data separation method step, similarity coefficient PJ satisfies above-mentioned functional expression (3.3), by average curved surface m _k(x, the two-dimentional error image h that y) separates _1k(x y) has reached the requirement of view data smooth change, then error image h _1k(x y) promptly is the steady data Layer image of the single yardstick information of being tried to achieve.

As mentioned above, for the two-dimensional image data matrix F (x, y) and error image h _{1 (k-1}) extreme point in (k 〉=2 o'clock) is distribution at random, it need be organized in an orderly manner by certain syntople, thereby form enveloping surface so that carry out the enveloping surface match.

According to the needs of enveloping surface match, can adopt the Delaunay triangulation, concrete steps are as follows:

1, determines divided region.Can form set M={P by the extreme point of image _i| i=1,2 ..., N} couples together the frontier point of image, forms line segment aggregate Q, asks the triangulation on the M on the zone that with Q is the border.

2, constitute optimum triangular shape according to optimizing criterion.For arbitrary line segment q among the Q _i, in M, having any and its formation optimum triangular shape at least, this just can be used as effective composition of subdivision result.Get q ₁=P ₁P ₂, in M, look for an extreme point P ₃, make P ₁P ₂P ₃There are not other extreme points in the triangle of 3 compositions, and require (| P ₁P ₂|-| P ₂P ₃|) ²+ (| P ₁P ₂|-| P ₁P ₃|) ²+ (| P ₂P ₂|-| P ₁P ₃|) ²Value is minimum, satisfies the P of above-mentioned condition ₃Point is unique confirmable, with P ₁P ₂Form local optimum triangle Δ ₁=Δ P ₁P ₂P ₃, this moment, divided region was reduced into G ₁=G-Δ P ₁P ₂P ₃

G ₁Again through same step definition G ₂=G ₁-Δ ₂, so continuing, divided region constantly dwindles, and is intact until subdivision.

Export the abutment points set of each extreme point, promptly other that constitute optimum triangular shape jointly with this extreme point have a few the set of composition, and the abutment points of all extreme points of image is gathered can constitute an adjacency list.

Definition syntople point_relation=(order, neighbour_point}).

Wherein first record is a discrete point numbering, and second record is the abutment points set of the corresponding discrete point of order institute behind the subdivision, in the set with the numbering replacement discrete point of discrete point.The all corresponding syntople of each discrete point, the union of all of its neighbor relation is formed adjacency list, is designated as relationship={point_relation_i|i=1, and 2 ..., N}.

When above-mentioned image extreme point carries out surface fitting, on the basis of the discrete point adjacency list that subdivision provided, adopt method of interpolation to guarantee second-order smooth on each bar adjacent side of triangle, the enveloping surface fitting effect that is obtained is good.

Adopt the operation steps of BB method of interpolation to be:

1, determines the adjoint point collection of arbitrary image extreme point P0.If extreme value point set M is last 1 P ₀(x ₀, y ₀) functional value be f ₀, can retrieve P according to the adjacency list of subdivision output ₀Adjoint point.If total m, form adjoint point collection NP={Q _i(x _i, y _i) | i=1,2 ..., m}.

2, on the P0 point definition quadratic polynomial function f (x, y)

f(x，y)＝c ₁(x-x ₀) ²+c ₂(y-y ₀) ²+c ₃(x-x ₀)(y-y ₀)

+c ₄(x-x ₀)+c ₅(y-y ₀)+f ₀

3, definition least square function

$g({\mathrm{<figure-callout\; id="1"\; label="c"\; filenames="C20041003633300211.png"\; state="\{\{state\}\}">c</figure-callout>}}_1,c_2,c_3,c_4,c_5)=\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}[f(x_i,y_i)-f_i{]}^2/d_i],$ Make it satisfy P ₀Point place's interpolation and about adjoint point Q _iThe least square approximation of ∈ NP.F wherein _iBe Q _iThe functional value that point is given, d _iBe Q _iTo P ₀Distance.

4, determine quadratic polynomial function f (x, coefficient y).By least square function character

$\frac{\&PartialD;g}{{\&PartialD;c}_i}=0,i=1,\mathrm{\&Lambda;},5$

Ask c _i, i=1 ..., 5.

5, the quadratic polynomial function f (x, y) determine after, obtain at P ₀One second derivative at some place.

6, structure BB interpolation polynomial.After the quadratic polynomial function that defines on each extreme point was determined, the value of any point P was determined by this leg-of-mutton three summits value and the position in triangle in the triangle.

If 3 P of non-colinear _i(x _i, y _i), i=1,2,3 form triangle, and (x, y) available triangle area coordinate system is expressed as its internal point P

P＝α ₁P ₁+α ₂P ₂+α ₃P ₃

0≤α wherein _i≤ 1, and $\underset{i=1}{\overset{3}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_i=1.$

Internal point P (x, value y) can be by interpolation polynomial: $f(x,y)=F\left(\mathrm{\&alpha;}\right)=\underset{\left|\mathrm{\&lambda;}\right|=n}{\mathrm{\&Sigma;}}{b}_{\mathrm{\&lambda;}}{B}_{\mathrm{\&lambda;}}^n\left(\mathrm{\&alpha;}\right)$ Obtain.

Wherein $B_{\mathrm{\&lambda;}}^n\left(\mathrm{\&alpha;}\right)=\frac{n!}{{\mathrm{\&lambda;}}_1!{\mathrm{\&lambda;}}_2!{\mathrm{\&lambda;}}_3!}{\mathrm{\&alpha;}}_1^{{\mathrm{\&lambda;}}_1}{\mathrm{\&alpha;}}_2^{{\mathrm{\&lambda;}}_2}{\mathrm{\&alpha;}}_3^{{\mathrm{\&lambda;}}_3}$ It is the Bernstein spline base function.

λ＝(λ ₁，λ ₂，λ ₃)，|λ|＝λ ₁+λ ₂+λ ₃。

When n=3, b is arranged ₃₀₀=f ₁, b ₀₃₀=f ₂, b ₀₀₃=f ₃, $b_{210}=f_1+\frac{1}{3}\&dtri;f{\left(P_1\right)}^T(P_2-P_1),$

$b_{201}={\mathrm{<figure-callout\; id="1"\; label="f"\; filenames="C20041003633300211.png"\; state="\{\{state\}\}">f</figure-callout>}}_1+\frac{1}{3}\&dtri;f{\left(P_1\right)}^T(P_3-P_1),$ $b_{120}=f_2+\frac{1}{3}\&dtri;f{\left(P_2\right)}^T(P_1-P_2),$

$b_{102}=f_3+\frac{1}{3}\&dtri;f{\left(P_3\right)}^T(P_1-P_3),$ $b_{021}=f_2+\frac{1}{3}\&dtri;f{\left(P_2\right)}^T(P_3-P_2),$

$b_{012}=f_3+\frac{1}{3}\&dtri;f{\left(P_3\right)}^T(P_2-P_3),$

$b_{111}={\mathrm{<figure-callout\; id="1"\; label="f"\; filenames="C20041003633300211.png"\; state="\{\{state\}\}">f</figure-callout>}}_1+\frac{1}{3}\&dtri;f{\left(P_1\right)}^T(P_2+P_3-2P_1)+\frac{1}{6}{(P_2-P_1)}^T{\&dtri;}^{2}f\left(P_1\right)(P_3-P_1).$

More than be with single figure layer h _1k(x is y) from original two dimensional image data matrix F (x, the basic performing step of separating in y).

Based on same principle, repeat above-mentioned single figure layer h _Mk(x, separating step y) n time can obtain the single figure layer that the n layer is separated from each other, and each layer result comprised the target information of single relatively yardstick.

(x isolates each single figure layer h in y) from original image F _Mk(x, y) after, should be no longer to contain isolated yardstick information in the residual image, the iteration stopping condition of its data separating method is:

Based on the error image h that judges by m layer average curved surface institute under true _1k, its corresponding similarity coefficient PJ all satisfies area boundary parameter alpha separately ₁, α ₂Determined function expression (3.3).

Or,

The average of coenvelope face that is made of the local maximum point and the lower envelope face that is made of local minizing point is zero, i.e. average curved surface m _{K (}X is zero y).

For ground floor, m=1.

H is arranged this moment _{1 (k-1)}(x, y)-m _1k(x, y)=h _Lk(x, y) (3.4)

Definition c ₁(x, y)=h _1k(x, y) (3.5)

This promptly is from isolated first figure layer information of original image, with c ₁(x, y) (x separates in y), and remaining part is designated as r from original SAR two dimensional image F ₁(x, y) then

F(x，y)-c ₁(x，y)＝r ₁(x，y) (3.6)

In fact, r ₁(x still might include a plurality of yardstick information in y).

With r ₁(x y) is the raw data that next figure layer data is separated, and repeats the single figure layer that the said method step can be isolated the m layer successively, and each layer result comprised the target information of single relatively yardstick.

As mentioned above, the 2-D data separation method of application electromagnetism echo has the following advantages and beneficial effect:

1, multiple dimensioned composite information treatment of picture effect has been optimized in the realization of two-dimentional EMD method greatly, for the information extraction of this type of image with analyze the disposal route that science is provided, make an amendment slightly and just may extend to three-dimensional;

2, pass through application experiment, a kind of really image processing method that has much development potentiality of two-dimentional as can be seen EMD method, in Military Early Warning, scouting, guidance, earth remote sensing observation in space is used, and use in national economy fields such as industry, medical science traffic bright prospects; Can be applicable to many aspects such as information separated, figure image intensifying, data tranquilization;

3, two-dimentional EMD method has adaptivity to the processing of different images, does not need manual intervention substantially, is a kind of posterior disposal route, and process can be accomplished object-oriented;

4, utilize two-dimentional EMD method that image is carried out yardstick and separate the characteristic of having utilized Riding Waves, very effective for the separation of images that overlapped information is arranged;

5, the border variation is all influential to the each iteration result in each layer, but influence is different.Pay attention to choosing the iterations of each layer in the two dimension EMD layering, reduced the influence of border the result.

Description of drawings

The present invention will be further described in conjunction with the following drawings.

Fig. 1 is the local data's extreme point distribution schematic diagram in the described two dimensional image;

Fig. 2-the 1st, inner neighborhood of a point synoptic diagram in a certain part in the two dimensional image;

Fig. 2-the 2nd, a certain border neighborhood of a point synoptic diagram in the two dimensional image;

Fig. 3 is behind the plane trigonometry subdivision one exemplary plot as a result;

Fig. 4 is the extremal surface exemplary plot that obtains by interpolation;

Fig. 5 is the schematic flow sheet of described two-dimensional image data separation method;

Fig. 6-the 1st, the sea SAR original image described in the embodiment 1;

Fig. 6-the 2nd, the tendency information that obtains at sea SAR separation of images to layer 5;

Fig. 6-the 3rd, the overlapped information of the h1-h5 that the sea SAR separation of images is come out.

Embodiment

Embodiment 1, and shown in Fig. 1 to Fig. 6-3, described two-dimensional image data separation method with multi-scale information is applied to survey the SAR treatment of picture of submarine topography information.

Shown in Fig. 6-1 is the sea SAR original two dimensional image of actual generation, wherein includes various features yardstick information, need therefrom progressively isolate the figure layer of the single yardstick of each self-separation such as ocean wave motion, ocean current field, wind field, mean field and circulation field.

At first isolate the hum pattern layer that embodies ocean wave motion from original SAR two dimensional image, as shown in Figure 5, the operation steps of described separation method is:

The first step is obtained above-mentioned SAR two-dimensional image data, promptly obtains the two-dimentional raw data in the image.

In second step, determine the extreme point in the above-mentioned SAR two dimensional image raw data.

Four drift angles in whole two-dimensional image data plane owing to consider the plane subdivision work whole plane of subdivision intactly that will guarantee next step, must force to be made as extreme point at the place, four summits of image.

Must cause data to have jiggly characteristic owing to survey the otherness of thing various piece character, therefore in the formed two-dimensional image data of electromagnetism echo plane, comprise a maximum point and a minimum point at least.If do not have the view data extreme point asked on the whole two dimensional surface, then by obtaining at least one maximum point and a minimum point after single order or several rank derivative operation.

Shown in Fig. 2-1, (i, j), next-door neighbour's point has 8 to the point of image inside.

Shown in Fig. 2-2, for the processing of data boundary, because it has only a half neighborhood that data are arranged, so can only seek extreme point in 1/2 neighborhood interval, next-door neighbour's point becomes 5.

As shown in Figure 1, be Local Extremum distribution schematic diagram in the electromagnetism echo SAR image.

In the 3rd step, extreme point is coupled together by certain syntople by the Delaunay triangulation.

After extreme point was found out, they were distribution at random in the plane, need organize them in an orderly manner by certain syntople, so that spatially carry out surface fitting.

Discrete point set M={Pi|i=1 on the given plane, 2 ..., N}, its outermost border point is coupled together, form a line segment aggregate Q={qi=PiPj|i ≠ j, 1≤i, j≤N, Pi, Pj ∈ M} requires to ask the triangulation on the M on the regional G that with Q is the border.

Select a kind of data structure of outstanding boundary sections.

If discrete point point=(order, x, y, z), wherein order is corresponding one by one with point, the numbering of expression point; (x, y z) are the locus of point;

Boundary sections edge=(order_1, order_2) boundary sections of expression divided region; Wherein, order_1, order_2 are the discrete point numberings, represent the line segment two-end-point.

(position_dist), position_dist is an orientation distance to definition phrase word=for index, order, is to be numbered the distance of the point of order to a certain point of fixity Po (being called anchor point), and index is the word index.

Each discrete point defines a corresponding word, and the opener of all word has constituted dictionary, is designated as dictionary={word _i| i=1,2 ..., N} is shown in table 2.2.1.

Described dictionary dictionary has roughly put down in writing the relative position relation of regional G point at random.By this dictionary dictionary,, can under the situation of 2 of known boundaries, optimized choice thirdly constitute optimum triangular shape, therefore can confine the position of optimum point basically.

Table 2.2.1

Output data adopts the adjacency list form of discrete point, is convenient to aftertreatment.Definition syntople point_relation=(order, neighbour_point}), wherein first record is a discrete point numbering, and second record is the abutment points set of the corresponding discrete point of order institute behind the subdivision, in the set with the numbering replacement discrete point of discrete point.

The all corresponding syntople of each discrete point, the union of all of its neighbor relation is formed adjacency list, is designated as relationship={point_relation_i|i=1, and 2 ..., N}, shown in table 2.2.2, it expresses the subdivision result.

Table 2.2.2

Export above-mentioned subdivision result and form the triangulation collection (x of rectangle _i ^j, y _i ^j, f _i ^j):

$\left(\begin{array}{ccccccccc}{x}_1^1& {\mathrm{<figure-callout\; id="1"\; label="y"\; filenames="C20041003633300211.png"\; state="\{\{state\}\}">y</figure-callout>}}_1^1& {\mathrm{<figure-callout\; id="1"\; label="f"\; filenames="C20041003633300211.png"\; state="\{\{state\}\}">f</figure-callout>}}_1^1& x_2^1& y_2^1& f_2^1& x_3^1& y_3^1& f_3^1\\ x_1^2& {\mathrm{<figure-callout\; id="1"\; label="y"\; filenames="C20041003633300211.png"\; state="\{\{state\}\}">y</figure-callout>}}_1^2& {\mathrm{<figure-callout\; id="1"\; label="f"\; filenames="C20041003633300211.png"\; state="\{\{state\}\}">f</figure-callout>}}_1^2& x_2^2& y_2^2& f_2^2& x_3^2& y_3^2& f_3^2\\ & .& .& .& .& .& .& & \\ x_1^m& y_1^m& f_1^m& x_2^m& y_2^m& f_2^m& x_3^m& y_3^m& f_3^m\end{array}\right)$

Above-mentioned triangulation is concentrated,

x _i ^jBe the horizontal ordinate on j leg-of-mutton i summit, j=1 wherein ..., m; M is leg-of-mutton total number; I=1,2,3.

y _i ^jBe the ordinate on j leg-of-mutton i summit, j=1 wherein ..., m; I=1,2,3.

f _i ^jBe (x _i ^j, y _i ^j) point of living in (i, given value j).

As shown in Figure 3, for arbitrary line segment qi among the Q, have any and its formation optimum triangular shape at least in M, this just can be used as effective composition of subdivision result.

Get q1=P ₁P ₂, optimize criterion according to certain, with the P among the M ₃Point is formed local optimum triangular shape Δ 1=Δ P ₁P ₂P ₃, this moment, divided region was reduced into G ₁=G-Δ P ₁P ₂P ₃

G ₁Again through same step definition G ₂=G ₁-Δ 2 so continues, and obtains following descending series: G ₀=G,

G ₁＝G ₀-Δ1，

G ₂，...，

Gn-1=Δ n until Gn=Gn-1-Δ n=Ф, has just finished the subdivision of regional G.

The 4th goes on foot, and tries to achieve the interpolation curved surface of average curved surface, and maximum point and the minimum point that exists in the curved surface carried out match and form its enveloping surface.

Described two-dimensional image data separation method, what extreme point was carried out that match adopts is the BB method of interpolation, is to carry out on the basis of the discrete point adjacency list that provides at subdivision.This method of interpolation can guarantee second-order smooth on each bar adjacent side of triangle, and is better to reach the surface fitting effect.

Set the P that discrete point set M goes up a bit ₀(x ₀, y ₀) functional value be f ₀, can retrieve P according to the adjacency list of subdivision output ₀Adjoint point.

Set above-mentioned P ₀The total m of adjoint point, then can form adjoint point collection NP={Q _i(x _i, y _i) | i=1,2 ..., m},

At P ₀Definition quadratic polynomial function on the point

f(x，y)＝c ₁(x-x ₀) ²+c ₂(y-y ₀) ²+c ₃(x-x ₀)(y-y ₀) (2.3.1)

+c ₄(x-x ₀)+c ₅(y-y ₀)+f ₀

Above-mentioned function (2.3.1) satisfies P ₀The interpolation at some place, but also satisfy about adjoint point Q _iThe least square approximation of ∈ NP.Promptly satisfy following least square function and reach minimum value.

$g({\mathrm{<figure-callout\; id="1"\; label="c"\; filenames="C20041003633300211.png"\; state="\{\{state\}\}">c</figure-callout>}}_1,c_2,c_3,c_4,c_5)=\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}[f(x_i,y_i)-f_i{]}^2/d_i]...(2.3.2)$

F wherein _iBe Q _iThe functional value that point is given, d _iBe Q _iTo P ₀Distance.Its reciprocal value is as Q _iPoint is to the weights of P point derivative value influence power.C then _i, i=1 ..., 5 satisfy following formula

$\frac{\&PartialD;g}{{\&PartialD;c}_i}=0,i=1,\mathrm{\&Lambda;},5---\left(\mathrm{2.3.3}\right)$

Above-mentioned functional expression (2.3.3) is launched arrangement is descended system of equations:

$\left(\begin{array}{ccccc}{D}_{11}& D_{12}& D_{13}& D_{14}& D_{15}\\ D_{21}& D_{22}& D_{23}& D_{24}& D_{25}\\ D_{31}& D_{32}& D_{33}& D_{34}& D_{35}\\ D_{41}& D_{42}& D_{43}& D_{44}& D_{45}\\ D_{51}& D_{52}& D_{53}& D_{54}& D_{55}\end{array}\right)\left(\begin{array}{c}{C}_1\\ C_2\\ C_3\\ C_4\\ C_5\end{array}\right)=\left(\begin{array}{c}{r}_1\\ r_2\\ r_3\\ r_4\\ r_5\end{array}\right)---\left(\mathrm{2.3.4}\right)$

The following formula brief note is: DC=R (2.3.5)

When above-mentioned matrix D was nonsingular matrix, this system of equations had unique solution.For guaranteeing that D is a nonsingular matrix, need expand the NP set sometimes, promptly comprise P ₀Non-direct adjoint point.They are at P ₀The approximate value of one second derivative at some place can specifically be expressed as follows:

$\&dtri;f\left(P\right){|}_{P=P_0}\&ap;{(\frac{\&PartialD;f}{\&PartialD;x},\frac{\&PartialD;f}{\&PartialD;y})}^T=\left(\begin{array}{c}{C}_4\\ C_5\end{array}\right)---\left(\mathrm{2.3.6}\right)$

${\&dtri;}^{2}f\left(P\right){|}_{P=P_0}{\&ap;\left(\begin{array}{cc}\frac{{\&PartialD;}^{2}f}{{\&PartialD;x}^2}& \frac{{\&PartialD;}^{2}f}{\&PartialD;x\&PartialD;y}\\ \frac{{\&PartialD;}^{2}f}{\&PartialD;x\&PartialD;y}& \frac{{\&PartialD;}^{2}f}{{\&PartialD;y}^2}\end{array}\right)}^T=\left(\begin{array}{cc}{2C}_1& C_3\\ C_3& {2C}_2\end{array}\right)---\left(\mathrm{2.3.7}\right)$

Following content is derived the BB interpolation polynomial.

Known triangle subdivision set T={tt _i| i=1,2 ..., M}, and vertex set V={ (x _i, y _i, f _i) ∈ R ³, construct Local C on each triangle ²(x y), satisfies the continuous interpolating function f of continuous and whole C1

f(x _i，y _i)＝f _i，i＝1，2，Λ，N (2.3.8)

Set 3 P of non-colinear _i(x _i, y _i), i=1,2,3 form its inside of triangle tt a bit.

(x, y) available triangle area coordinate system is expressed as P

P＝α ₁P ₁+α ₂P ₂+α ₃P ₃ (2.3.9)

0≤α wherein _i≤ 1, and $\underset{i=1}{\overset{3}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_i=1.$

Set $f(x,y)=\underset{l=0}{\overset{n}{\mathrm{\&Sigma;}}}\underset{i+j=l}{\mathrm{\&Sigma;}}{a}_{\mathrm{ij}}{x}^i{y}^i---\left(\mathrm{2.3.10}\right)$

Wherein equation (2.3.9) deformable is

$\left(\begin{array}{c}x\\ y\end{array}\right)={\mathrm{\&alpha;}}_i\left(\begin{array}{c}{x}_1\\ {\mathrm{<figure-callout\; id="1"\; label="y"\; filenames="C20041003633300211.png"\; state="\{\{state\}\}">y</figure-callout>}}_1\end{array}\right)+{\mathrm{\&alpha;}}_2\left(\begin{array}{c}{x}_2\\ y_2\end{array}\right)+{\mathrm{\&alpha;}}_3\left(\begin{array}{c}{x}_3\\ y_3\end{array}\right)---\left(\mathrm{2.3.11}\right)$

With formula (2.3.11) substitution formula (2.3.10), draw through conversion:

$f(x,y)=F\left(\mathrm{\&alpha;}\right)=\underset{\left|\mathrm{\&lambda;}\right|=n}{\mathrm{\&Sigma;}}{b}_{\mathrm{\&lambda;}}{B}_{\mathrm{\&lambda;}}^n\left(\mathrm{\&alpha;}\right)---\left(\mathrm{2.3.12}\right)$

Wherein, $B_{\mathrm{\&lambda;}}^n\left(\mathrm{\&alpha;}\right)=\frac{n!}{{\mathrm{\&lambda;}}_1!{\mathrm{\&lambda;}}_2!{\mathrm{\&lambda;}}_3!}{\mathrm{\&alpha;}}_1^{{\mathrm{\&lambda;}}_1}{\mathrm{\&alpha;}}_2^{{\mathrm{\&lambda;}}_2}{\mathrm{\&alpha;}}_3^{{\mathrm{\&lambda;}}_3}$ It is the Bernstein spline base function.

λ＝(λ ₁，λ ₂，λ ₂)，|λ|＝λ ₁+λ ₂+λ ₃。

When n=3, b is arranged ₃₀₀=f ₁, b ₀₃₀=f ₂, b ₀₀₃=f ₃, $b_{210}={\mathrm{<figure-callout\; id="1"\; label="f"\; filenames="C20041003633300211.png"\; state="\{\{state\}\}">f</figure-callout>}}_1+\frac{1}{3}\&dtri;f{\left(P_1\right)}^T(P_2-P_1),$

$b_{201}={\mathrm{<figure-callout\; id="1"\; label="f"\; filenames="C20041003633300211.png"\; state="\{\{state\}\}">f</figure-callout>}}_1+\frac{1}{3}\&dtri;f{\left(P_1\right)}^T(P_3-P_1),$ $b_{120}=f_2+\frac{1}{3}\&dtri;f{\left(P_2\right)}^T(P_1-P_2),$

$b_{102}=f_3+\frac{1}{3}\&dtri;f{\left(P_3\right)}^T(P_1-P_3),$ $b_{021}=f_2+\frac{1}{3}\&dtri;f{\left(P_2\right)}^T(P_3-P_2),$

$b_{012}=f_3+\frac{1}{3}\&dtri;f{\left(P_3\right)}^T(P_2-P_3),$

$b_{111}={\mathrm{<figure-callout\; id="1"\; label="f"\; filenames="C20041003633300211.png"\; state="\{\{state\}\}">f</figure-callout>}}_1+\frac{1}{3}\&dtri;f{\left(P_1\right)}^T(P_2+P_3-2P_1)+\frac{1}{6}{(P_2-P_1)}^T{\&dtri;}^{2}f\left(P_1\right)(P_3-P_1).$

In sum, according to the enveloping surface that the determined function of equation (2.3.12) is derived, just try to achieve the interpolation curved surface of average curved surface.

In the 5th step, from original SAR two dimensional image, isolate the target information figure layer that embodies ocean wave motion.

Through above-mentioned first to fourth step, each Local Extremum has been carried out the Delaunay triangulation, and each subdivision triangle is carried out interpolation to guarantee second-order smooth on the border, finally obtains F (x, the coenvelope face X that is determined by maximum point y) _Max(x, y), and by the definite lower envelope face X of minimum point _Min(x, y).

By coenvelope face X _Max(x is y) with lower envelope face X _Min(x y) obtains average curved surface, i.e. a m ₁(x, y)=[X _Max(x, y)+X _Min(x, y)]/2.

Average curved surface m ₁(x, y) difference with original digital image data is designated as h ₁(x, y), that is: h ₁(x, y)=F (x, y)-m ₁(x, y) (3.1)

Because h ₁(x, waveform y) still includes the phenomenon of compound fluctuation, and the local mean value curved surface is also also non-vanishing, and therefore repeating above-mentioned first to fourth step obtains next average curved surface once more, and h ₁(x is y) as the raw data of next average curved surface.

Just except when for the first time trying to achieve the average curved surface, be that (x, y) as beyond the raw data, later all is that (x y) is raw data for poor h with the last time when repeatedly trying to achieve the average curved surface with two-dimensional SAR image data F.

Then secondary difference is designated as,

h ₁(x，y)-m ₁₁(x，y)＝h ₁₁(x，y) (3.2)

Repeat above-mentioned steps k time, up to trying to achieve average curved surface m _k(x is y) with the determined poor h of original image _1k, try to achieve average curved surface m with the last time _K-1(x, y) fixed poor h _{1 (k-1)}The standard-required that meets similarity coefficient PJ.

Wherein, similarity coefficient PJ satisfies following functional expression.

$\mathrm{PJ}=\underset{x=1}{\overset{W}{\mathrm{\&Sigma;}}}\underset{y=1}{\overset{H}{\mathrm{\&Sigma;}}}\left[\frac{{|(h_{m(k-1}(x,y)-h_{\mathrm{mk}}(x,y)|}^2}{h_{m(k-1)}^2(x,y)}\right]\&Element;[{\mathrm{\&alpha;}}_1,{\mathrm{\&alpha;}}_2]---\left(3.3\right)$

In above-mentioned functional expression,

M is the number of plies (m=1 in fact at this moment) when pre-treatment;

K is the number of times of this layer re-treatment;

W and H are respectively picture traverse and height;

And α ₁And α ₂Be the steady data Layer area boundary parameter of given embodiment ocean wave motion information, α ₁And α ₂Be constant, α ₁＜α ₂

If similarity coefficient PJ satisfies above-mentioned functional expression (3.3), the average curved surface m that then tries to achieve _k(x is y) with the last m that tries to achieve _K-1(x y) compares definite similarity coefficient PJ, is to be in the area boundary parameter alpha ₁And α ₂Between, the h that is tried to achieve then is described _1k(x y) has reached the requirement of view data smooth change, h _1k(x, the steady data Layer image of the embodiment ocean wave motion information of y) being tried to achieve just.

Foregoing is based on the data separating process of a dimensional properties of ocean wave motion.

Based on design of same method and principle, as shown in Figure 5, isolate the data of ground floor ocean wave motion yardstick information after, can isolate the single figure layer data of each self-separation that comprises ocean current field, wind field, mean field and circulation field step by step.This moment is with the area boundary parameter alpha ₁, α ₂The eigenwert that is set at ocean current field, wind field, mean field or circulation field respectively gets final product.

Try to achieve average curved surface m _k(x, iteration stopping condition y):

Generally be the poor h when trying to achieve m layer average curved surface based on judgement _1k, similarity coefficient PJ all satisfies area boundary parameter alpha separately ₁, α ₂The function expression of determining (3.3).

Or the average of coenvelope face that is made of the local maximum point and the lower envelope face that is made of local minizing point is zero.

For ground floor, m=1.

H is arranged this moment _{1 (k-1)}(x, y)-m _1k(x, y)=h _1k(x, y) (3.4)

Definition c ₁(x, y)=h _1k(x, y) (3.5)

This promptly is from isolated first figure layer information of original image, with c ₁(x, y) (x separates in y), and remaining part is designated as r from original SAR two dimensional image F ₁(x, y) then

F(x，y)-c ₁(x，y)＝r ₁(x，y) (3.6)

In fact, r ₁(x still includes yardstick information such as ocean current field, wind field, mean field and circulation field in y).

With r ₁(x is the raw data that next figure layer data is separated y), repeats the above-mentioned method step that comprises the steady data Layer of ocean wave motion information, can isolate the data plot layer of the single yardstick information of ocean current field, wind field, mean field or circulation field successively.

In theory, the two-dimensional image data that includes a plurality of yardstick information separates and can operate n time, represents with following iterative relation equation:

r ₁(x，y)-c ₂(x，y)＝r ₂(x，y)，…，r _n-1(x，y)-c _n(x，y)＝r _n(x，y) (3.7)

The termination condition of above-mentioned n layer average curved surface decomposable process is:

Work as r _n(x, y) or c _n(x, y) determined similarity coefficient PJ all satisfies yardstick information area limit parameter (α separately ₁, α ₂), as function expression (3.3);

Or, work as r _n(x, y) become a smooth change figure layer and can't be again when wherein obtaining the average curved surface till, then final remainder just will be antiderivative trend term.

Consolidated equation formula (3.6) and (3.7) obtain following function expression:

$F(x,y)=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{h}_i(x,y)+r_n(x,y)---\left(3.8\right)$

Wherein, h _i(x y) is the different scale figure layer that obtains after separating according to yardstick, r _n(x y) is the final trend term that obtains.

As the above embodiments content, use two-dimensional image data separation method treatment S AR view data, can go out the single figure layer and the final trend term of each self-separation of ocean wave motion, ocean current field, wind field, mean field and circulation field according to the different scale information separated, for the research of each yardstick phenomenon provides data, also particularly the research of image analysis of spectrum and Study on Target Recognition provide a kind of active data source preprocess method for the ocean remote sensing applied research.

As Fig. 6-the 2nd, the sea SAR original image utilizes this method to handle and is separated to the tendency information that obtains behind the layer 5, can find out the distribution situation on ocean and land from image significantly.

As Fig. 6-the 3rd, the image information that obtains after the h1-h5 stack that the sea SAR separation of images is come out.

Claims (5)

1, a kind of two-dimensional image data separation method with multi-scale information, it is characterized in that: the execution in step of described separation method is, at two-dimensional image data matrix F (x, y) determine the extreme point that each is local, and couple together by certain syntople and to form the subdivision triangle, the interpolation method of employing second-order smooth on each bar adjacent side of triangle carries out interpolation, carry out surface fitting for maximum point in the image and minimum point, thereby form coenvelope face Xmax (x, y) and lower envelope face Xmin (x, y), obtain an average curved surface m by upper and lower envelope surface ₁(x, y);

(x gets rid of average curved surface m in y) from original image F _k(x y) obtains error image h _1k, work as k=1, promptly by average curved surface m ₁(x y) obtains error image h ₁₁

With error image h ₁₁Be the original image of data separating next time, repeat above-mentioned detachment process and obtain average curved surface m ₂(x is y) with error image h ₁₂, obtain average curved surface m until repeating above-mentioned separation method k time _k(x is y) with error image h _1k

Wherein, the k secondary data is separated the error image h that obtains _1kBe an above error image h _{1 (k-1)}For original image, wherein k is constant and k 〉=2.

As difference image h _1kThe average curved surface m that obtains with the last time _K-1(x, y) determined error image h _{1 (k-1)}, when satisfying following similarity coefficient PJ function expression, above-mentioned data separating process end:

$\mathrm{PJ}=\underset{x=1}{\overset{W}{\mathrm{\&Sigma;}}}\underset{y=1}{\overset{H}{\mathrm{\&Sigma;}}}\left[\frac{{|(h_{m(k-1)}(x,y)-h_{\mathrm{mk}}(x,y)|}^2}{h_{m(k-1)}^2(x,y)}\right]\&Element;[{\mathrm{\&alpha;}}_1,{\mathrm{\&alpha;}}_2]---\left(3.3\right)$

Wherein, m is the number of plies when pre-treatment;

K is the number of times of this layer re-treatment, and k 〉=2;

W and H are respectively picture traverse and height;

And α ₁And α ₂Be carry out before the data separating the steady data Layer area boundary parameter of given single yardstick information, α ₁And α ₂Be constant, and α ₁＜α ₂

Error image h _Mk(x y) promptly is the mask data image with single yardstick information.

2, the two-dimensional image data separation method with multi-scale information according to claim 1 is characterized in that: (x isolates single figure layer h in y) from original image F _Mk(x, y) after, with residual image r _m(x y) is the raw data that next figure layer data is separated;

Repeat above-mentioned steps m time, can isolate the single figure layer of m layer successively;

The iteration stopping condition of described data separating method is,

Based on judging by the determined error image h of m layer average curved surface _1k, its corresponding similarity coefficient PJ all satisfies area boundary parameter alpha separately ₁, α ₂Determined function expression 3.3;

Or the average of coenvelope face that is made of the local maximum point and the lower envelope face that is made of local minizing point is zero, i.e. average curved surface m _k(x is zero y).

3, the two-dimensional image data separation method with multi-scale information according to claim 2 is characterized in that: for the two-dimensional image data matrix F (x, y) and error image h _{1 (k-1)}In extreme point, adopt triangulation that extreme point is organized in an orderly manner by certain syntople, so that carry out the enveloping surface match.

4, the two-dimensional image data separation method with multi-scale information according to claim 3, it is characterized in that: when carrying out the enveloping surface match, on the basis of the discrete point adjacency list that subdivision provided, adopt the method for interpolation of second-order smooth on each bar adjacent side of triangle.

5, according to claim 3 or 4 described two-dimensional image data separation methods with multi-scale information, it is characterized in that: (x is y) with error image h in definite above-mentioned two-dimensional image data matrix F _{1 (k-1)}In extreme point the time, four drift angles in the whole two-dimensional image data plane are made as extreme point.

Images