CN101676745B - Method for separating gravity and magnetic field in DCT domain - Google Patents

Abstract

The invention belongs to the field of processing of gravity and magnetic data, in particular discloses a method for separating a gravity and magnetic field in a DCT domain, which comprises the following steps: carrying out gridding treatment on two-dimensional discrete gravity and magnetic data obtained by field measurement, and forming regular grid data; converting the grid data into the DCT domain from a space domain, and obtaining a frequency spectrum of original data; setting n different threshold values En, and counting an area An surrounded by all sets, the energy spectrum values of which are larger than En; drawing a logEn-logAn bilogarithmic graph in a rectangular Cartesian coordinate system, using a least square method to fit a plurality of straight-line segments with different slopes, and determining an energy spectrum boundary value of a filter according to the energy spectrum value corresponding to the intersection point of the straight-line segments; carrying out filtering treatment on the gravity and magnetic field according to the determined energy spectrum boundary value; and carrying out the inverse discrete cosine transform to the filtering result and converting the filtering result into the space domain from the DCT domain to acquire the separating result of the original data. The method of the invention has stronger decorrelation ability, fully considers the spatial self-similarity or multi-fractal characteristics of the gravity and magnetic field, and can effectively separate the gravity and magnetic field from a local field.

Description

Method for separating gravity and magnetic field in a kind of DCT territory

Technical field

The invention belongs to heavy magnetic data process field, be specifically related to the method for separating gravity and magnetic field in a kind of DCT territory.

Background technology

Zone be one of important content of heavy magnetic data work of treatment separating of local field.The 1950's, the proposition of discrete fourier transition (DFT), the particularly appearance of nineteen sixty-five rapid fourier change algorithm (FFT), heavy magnetic data is handled breaks away from from the spatial domain computing of complexity, also make filtering method in the DFT territory become to attach most importance to main method that magnetic field separates is as high pass commonly used, low-pass filtering etc.These filtering methods carry out based on the size of frequency often, but in fact similar geologic body also can have different frequencies, on the contrary, therefore geophysical field with same frequency also may not reflect similar geological phenomenon, traditional filtering method geophysical field of being difficult to distinguish similar geologic body or having identical geological Significance.1975, fractal theory that Mandelbrot proposes and the multifractal theory that grew up afterwards were for the statistical property of studying frequency spectrum in frequency field provides theoretical foundation.A large amount of examples show that many geological processes comprise that mineralization process all can cause the field-multifractal field with self-similarity or statistical self-similarity.Therefore, the ripple filter is handled except that the frequency difference of considering the field, also should pay attention to the space self-similarity of geophysical field.Based on above background, Cheng Qiuming proposed the multifractal filtering method (" power spectrum-area " is called for short the S-A method) in the DFT frequency domain in 1985, and the decomposition or the further feature that are used for spatial model (as unusual and background) extract.This method has been used to the decomposition of geochemical anomaly and remote sensing images etc. at present.

In fact, in heavy magnetic measurement data, the field, zone is always overlapped with the frequency spectrum of local field, the correlativity of various degrees between them.How reducing effectively or remove correlativity between the heavy magnetic field of different scale, thereby carry out accurately and effectively separating of field, zone and local field, is the important topic that the geophysics worker faces.The discrete cosine transform (DCT) that three India descendants such as Ahmed American scholar proposed in 1974 has obtained at aspects such as voice, picture coding and data compressions using widely with its good performance.DCT is the same with DFT, all belongs to the sinusoidal class conversion in the orthogonal transformation, and the condition of its existence is identical with the Fourier integral condition of convergence, and has the character similar to DFT in some aspects.Different is, DCT has the ability of the decorrelation stronger than DFT, reason is that DCT has very strong " concentration of energy " characteristic: the energy of most of natural signs (comprising sound and image etc.) all concentrates on the low frequency part in DCT territory, and when signal had statistical property near Markov process (Markovprocesses), the decorrelation of DCT approached to have the performance of the Karhunen-Loeve transformation (Karhunen-Loeve) of optimum decorrelation.DFT requires signal only can have limited discontinuous point and extreme point in one-period or in the limit of integration, and definitely can amass (satisfying the Di Lihelai condition).Because the existence in implicit cycle can cause the discontinuous of edge among the DFT, be subjected to the influence of Gibbs' effect easily and make the result of limit portion produce bigger distortion.In addition, for real continuous signal, DCT can avoid a large amount of complex operations, and for DFT, operation efficiency is higher and be easy to programming and realize.

The application of DCT both at home and abroad at present mainly is confined to fields such as geological data processing and compression of images.2005, Zhang Fengxu etc. are applied to DCT in the heavy magnetic data processing first, study the conventional treatment methods such as upward continuation, differentiate of gravity field in the DCT territory, but do not related to space self-similarity or the multi-fractal features and the associated method for separating gravity and magnetic field in heavy magnetic field.

Summary of the invention

The object of the present invention is to provide the method for separating gravity and magnetic field in a kind of DCT territory, can separate heavy magnetic area field and local field effectively.

Realize the technical scheme of the object of the invention, the method for separating gravity and magnetic field in a kind of DCT territory may further comprise the steps:

(1) the heavy magnetic measurement data of the two-dimensional discrete that will adopt heavy magnetic measurement instrument to measure are in the open air carried out gridding processing, the grid data of formation rule;

(2) utilize forward discrete cosine transform (DCT), the grid data that forms in the step (1) is changed in the DCT territory from transform of spatial domain, obtain the frequency spectrum data of raw data in the DCT territory; Described step

Formula used when carrying out forward discrete cosine transform (2) is shown in following formula (1):

$C(u,v)=a\left(u\right)a\left(v\right)\underset{m=0}{\overset{M-1}{\mathrm{\Σ}}}\underset{n=0}{\overset{N-1}{\mathrm{\Σ}}}f(m,n)\cos \frac{\mathrm{\π}(2m+1)u}{2M}\cos \frac{\mathrm{\π}(2n+1)v}{2N}$ Formula (1)

In the formula (1), f (m, n): m=0,1 ..., M-1; N=0,1 ..., N-1} is the heavy magnetic measurement data of the two dimension in the spatial domain, M, N are respectively the line number and the columns of data, C (u, v): u=0,1 ..., M-1; V=0,1 ..., the attach most importance to forward discrete cosine transform value of magnetic measurement data of N-1}, conversion coefficient a (u) and a (v) be defined as follows shown in the formula (2) of face:

$a\left(u\right)=\left\{\begin{array}{cc}1/\sqrt{M},& u=0\\ \sqrt{2/M},& 1\≤u\≤M-1\end{array}\right.,$ $a\left(v\right)=\left\{\begin{array}{cc}1/\sqrt{N},& v=0\\ \sqrt{2/N},& 1\≤v\≤N-1\end{array}\right.$ Formula (2)

(3) in the DCT territory, set n different threshold value E _n, these threshold values are between the maximal value and minimum value of power spectrum; For each threshold value E _n, the area A that statistics power spectrum value is surrounded greater than its all set _nThe computing formula of power spectrum is shown in following formula (3) in the middle DCT territory of described step (3):

E=|C (u, v) | ²Formula (3)

According to the maximal value and the minimum value of power spectrum, determine n different threshold value E _n, n is the integer between 100～1000, the power spectrum value is greater than E in the statistics DCT territory _nThe grid number, the product of grid number and grid area promptly is that the power spectrum value is greater than E _nAll area A of being surrounded of set _n

(4) in Descartes's rectangular coordinate system, respectively with LogE _nAnd LogA _nBe horizontal ordinate and ordinate, draw bilogarithmic graph;

(5) on the bilogarithmic graph of drawing in step (4), utilize least square fitting to go out many straight-line segments that slope is different, the power spectrum value of straight-line segment intersection point place correspondence is used for determining the power spectrum cut off value of different wave filters;

(6) according to the power spectrum cut off value of determining in the step (5), Filtering Processing is carried out in counterweight magnetic field in the DCT territory;

(7) the filtering result to obtaining in the step (6) carries out inverse discrete cosine transform (IDCT), and it is transformed in the spatial domain from the DCT territory, obtains the separating resulting of raw data; Formula used when carrying out inverse discrete cosine transform in the described step (7) is shown in following formula (4):

$f(m,n)=a\left(u\right)a\left(v\right)\underset{m=0}{\overset{M-1}{\mathrm{\Σ}}}\underset{n=0}{\overset{N-1}{\mathrm{\Σ}}}C(u,v)\cos \frac{\mathrm{\π}(2m+1)u}{2M}\cos \frac{\mathrm{\π}(2n+1)v}{2N}$ Formula (4)

Identical in every implication and formula (1), (2) in the formula (4).

The logarithm that uses in the described step (4) is the end with natural number, in Descartes's rectangular coordinate system, draws LogE _n-LogA _nBilogarithmic graph.

The principle of filtering is in the described step (6): according to the power spectrum cut off value of the field of determining, zone with local field, when separating local field, with the C (u of power spectrum value in the DCT territory greater than all set of power spectrum cut off value, v) the value tax is zero, otherwise, when the separated region field, (u, v) to compose be zero to value less than the C of all set of power spectrum cut off value with power spectrum value in the DCT territory.

Effect of the present invention is: method of the present invention has very strong decorrelation ability, has taken into full account the space self-similarity or the multi-fractal features in heavy magnetic field simultaneously, and separate with local field counterweight magnetic area field effectively.

Fig. 1 is the method for separating gravity and magnetic field process flow diagram in a kind of DCT provided by the present invention territory.

Fig. 2 is grouan area, a Jiangxi boat magnetic Δ T abnormal data equal-value map.

Description of drawings

Fig. 3 is the power spectrum-area bilogarithmic graph of boat magnetic Δ T data in the DCT territory.

Boat magnetic area field, grouan area, the Jiangxi equal-value map that Fig. 4 obtains for adopting the present invention to separate.

Grouan area, the Jiangxi boat magnetic local field equal-value map that Fig. 5 obtains for adopting the present invention to separate.

Below in conjunction with the drawings and specific embodiments the method for separating gravity and magnetic field in a kind of DCT provided by the present invention territory is described in further detail.

As shown in Figure 1, the method for separating gravity and magnetic field in a kind of DCT territory may further comprise the steps:

Embodiment

(1) the heavy magnetic measurement data of two-dimensional discrete that will adopt conventional heavy magnetic measurement instrument to measure are in the open air carried out gridding processing, the grid data of formation rule;

(2) utilize forward discrete cosine transform (DCT), the grid data that forms in the step (1) is changed in the DCT territory from transform of spatial domain, obtain the frequency spectrum data of raw data in the DCT territory;

Formula used when carrying out forward discrete cosine transform in the step (2) is shown in following formula (1):

$C(u,v)=a\left(u\right)a\left(v\right)\underset{m=0}{\overset{M-1}{\mathrm{\Σ}}}\underset{n=0}{\overset{N-1}{\mathrm{\Σ}}}f(m,n)\cos \frac{\mathrm{\π}(2m+1)u}{2M}\cos \frac{\mathrm{\π}(2n+1)v}{2N}$ Formula (1)

In the formula (1), f (m, n): m=0,1 ..., M-1; N=0,1 ..., N-1} is the heavy magnetic measurement data in the spatial domain, M, N are respectively the line number and the columns of data, C (u, v): u=0,1 ..., M-1; V=0,1 ..., N-1} is the forward discrete cosine transform value, conversion coefficient a (u) and a (v) be defined as follows shown in the formula (2) of face:

$a\left(u\right)=\left\{\begin{array}{cc}1/\sqrt{M},& u=0\\ \sqrt{2/M},& 1\≤u\≤M-1\end{array}\right.,$ $a\left(v\right)=\left\{\begin{array}{cc}1/\sqrt{N},& v=0\\ \sqrt{2/N},& 1\≤v\≤N-1\end{array}\right.$ Formula (2)

(3) in the DCT territory, set n different threshold value E _n, these threshold values are between the maximal value and minimum value of power spectrum; For each threshold value E _n, the area A that statistics power spectrum value is surrounded greater than its all set _n

The computing formula of power spectrum is shown in following formula (3) in the middle DCT territory of step (3):

E=|C (u, v) | ²Formula (3)

According to the maximal value and the minimum value of power spectrum, determine n different threshold value E _n, n is the integer between 100～1000, the power spectrum value is greater than E in the statistics DCT territory _nThe grid number, the product of grid number and grid area promptly is that the power spectrum value is greater than E _nAll area A of being surrounded of set _n

(4) in Descartes's rectangular coordinate system, respectively with LogE _nAnd LogA _nBe horizontal ordinate and ordinate, draw bilogarithmic graph;

The logarithm that uses in the step (4) is the end with natural number, in Descartes's rectangular coordinate system, draws LogE _n-LogA _nBilogarithmic graph.

(5) on the bilogarithmic graph of drawing in step (4), utilize least square fitting to go out many straight-line segments that slope is different, the power spectrum value of straight-line segment intersection point place correspondence is used for determining the power spectrum cut off value of different wave filters;

(6) according to the power spectrum cut off value of determining in the step (5), Filtering Processing is carried out in counterweight magnetic field in the DCT territory;

The principle of filtering is in the step (6): according to the power spectrum cut off value of the field of determining, zone with local field, when separating local field, with the C (u of power spectrum value in the DCT territory greater than all set of power spectrum cut off value, v) the value tax is zero, otherwise, when the separated region field, (u, v) to compose be zero to value less than the C of all set of power spectrum cut off value with power spectrum value in the DCT territory.

(7) the filtering result to obtaining in the step (6) carries out inverse discrete cosine transform (IDCT), and it is transformed in the spatial domain from the DCT territory, obtains the separating resulting of raw data.

Formula used when carrying out inverse discrete cosine transform in the step (7) is shown in following formula (4):

$f(m,n)=a\left(u\right)a\left(v\right)\underset{m=0}{\overset{M-1}{\mathrm{\Σ}}}\underset{n=0}{\overset{N-1}{\mathrm{\Σ}}}C(u,v)\cos \frac{\mathrm{\π}(2m+1)u}{2M}\cos \frac{\mathrm{\π}(2n+1)v}{2N}$ Formula (4)

Identical in every implication and formula (1), (2) in the formula (4).

As shown in Figure 2, be example with the boat magnetic Δ T abnormal data in grouan area, Jiangxi, the boat magnetic area field and the local field in this district separated.The present invention handles boat magnetic data according to described 7 steps of summary of the invention, that is: when indoor processing

(1) the boat magnetic Δ T abnormal data of the two-dimensional discrete that field survey is obtained carries out the gridding processing, forms 561 * 505 grid data, and the some line-spacing is 0.1km;

(2) utilize forward discrete cosine transform (DCT), the grid data that forms in the step (1) is changed in the DCT territory from transform of spatial domain, obtain the frequency spectrum data of raw data in the DCT territory;

(3) in the DCT territory, set n different threshold value E _n, these threshold values are between the maximal value and minimum value of power spectrum; For each threshold value E _n, the area A that statistics power spectrum value is surrounded greater than its all set _n

(4) in Descartes's rectangular coordinate system, respectively with LogE _nAnd LogA _nBe horizontal ordinate and ordinate, draw bilogarithmic graph, as shown in Figure 3;

(5) on the bilogarithmic graph of drawing in step (4), utilize least square fitting to go out many straight-line segments that slope is different, the power spectrum value of straight-line segment intersection point place correspondence is used for determining the power spectrum cut off value of different wave filters;

(6), in the DCT territory, boat magnetic Δ T abnormal data is carried out the Filtering Processing of zone and local field respectively according to determined power spectrum cut off value in the step (5);

(7) respectively to the zone that obtains in the step (6) and the filtering result of local field, carry out inverse discrete cosine transform (IDCT), it is transformed in the spatial domain from the DCT territory, magnetic area field (as shown in Figure 4) and local field (as shown in Figure 5) obtain navigating.

Claims (3)

1. the method for separating gravity and magnetic field in the DCT territory is characterized in that: may further comprise the steps:

(1) the heavy magnetic measurement data of the two-dimensional discrete that will adopt heavy magnetic measurement instrument to measure are in the open air carried out gridding processing, the grid data of formation rule;

(2) utilize forward discrete cosine transform (DCT), the grid data that forms in the step (1) is changed in the DCT territory from transform of spatial domain, obtain the frequency spectrum data of raw data in the DCT territory; Formula used when carrying out forward discrete cosine transform in the described step (2) is shown in following formula (1):

Formula (1)

In the formula (1), F (m, n): m=0,1 ..., M-1; N=0,1 ..., N-1Be the heavy magnetic measurement data of two dimension in the spatial domain, M, N are respectively the line number and the columns of data, C (u, v): u=0,1 ..., M-1; V=0,1 ..., N-1The forward discrete cosine transform value of the magnetic measurement data of attaching most importance to, conversion coefficient A (u)With A (v)Be defined as follows shown in the formula (2) of face:

Formula (2)

(3) in the DCT territory, set n different threshold value E _n, these threshold values are between the maximal value and minimum value of power spectrum, for each threshold value E _n, the area A that statistics power spectrum value is surrounded greater than its all set _nThe computing formula of power spectrum is shown in following formula (3) in the middle DCT territory of described step (3):

E=|C (u, v) | ² Formula (3)

According to the maximal value and the minimum value of power spectrum, determine n different threshold value E _n, n is the integer between 100～1000; The power spectrum value is greater than E in the statistics DCT territory _nThe grid number, the product of grid number and grid area promptly is that the power spectrum value is greater than E _nAll area A of being surrounded of set _n

(4) in Descartes's rectangular coordinate system, respectively with LogE _nAnd LogA _nBe horizontal ordinate and ordinate, draw bilogarithmic graph;

(5) on the bilogarithmic graph of drawing in step (4), utilize least square fitting to go out many straight-line segments that slope is different, the power spectrum value of straight-line segment intersection point place correspondence is used for determining the power spectrum cut off value of different wave filters;

(6) according to the power spectrum cut off value of determining in the step (5), Filtering Processing is carried out in counterweight magnetic field in the DCT territory;

(7) the filtering result to obtaining in the step (6) carries out inverse discrete cosine transform (IDCT), and it is transformed in the spatial domain from the DCT territory, obtains the separating resulting of raw data; Formula used when carrying out inverse discrete cosine transform in the described step (7) is shown in following formula (4):

Formula (4)

Identical in every implication and formula (1), (2) in the formula (4).

2. the method for separating gravity and magnetic field in a kind of DCT according to claim 1 territory is characterized in that: the logarithm that uses in the described step (4) is the end with natural number, in Descartes's rectangular coordinate system, draws LogE _n-LogA _nBilogarithmic graph.

3. the method for separating gravity and magnetic field in a kind of DCT according to claim 1 territory, it is characterized in that: the principle of filtering is in the described step (6): according to the power spectrum cut off value of the field of determining, zone with local field, when separating local field, with power spectrum value in the DCT territory greater than all set of power spectrum cut off value C (u, v)It is zero that value is composed, otherwise, when the separated region field, with power spectrum value in the DCT territory less than all set of power spectrum cut off value C (u, v)It is zero that value is composed.

Images