CN107632964A - A kind of plane GEOMAGNETIC FIELD downward continuation recurrence cosine transform method - Google Patents
A kind of plane GEOMAGNETIC FIELD downward continuation recurrence cosine transform method Download PDFInfo
- Publication number
- CN107632964A CN107632964A CN201710795560.4A CN201710795560A CN107632964A CN 107632964 A CN107632964 A CN 107632964A CN 201710795560 A CN201710795560 A CN 201710795560A CN 107632964 A CN107632964 A CN 107632964A
- Authority
- CN
- China
- Prior art keywords
- continuation
- formula
- geomagnetic field
- plane
- cosine transform
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Landscapes
- Measuring Magnetic Variables (AREA)
- Complex Calculations (AREA)
Abstract
The invention discloses a kind of plane GEOMAGNETIC FIELD downward continuation recurrence cosine transform method, belong to earth magnetism potential field data processing field.Concretely comprise the following steps：The pretreatment of plane GEOMAGNETIC FIELD measurement data；The gridding processing of plane GEOMAGNETIC FIELD measurement data；Utilize the cosine transform wave spectrum of recurrence cosine transform Calculation Plane GEOMAGNETIC FIELD；One a small amount of ε is set according to the requirement of downward continuation precision, chooses damping factor α value；The wave spectrum of continuation plane GEOMAGNETIC FIELD field is directly calculated by the compensating factor of inspection surface cosine transform wave spectrum and plane GEOMAGNETIC FIELD downward continuation；Wave spectrum is entered every trade to row to one-dimensional recurrence anticosine become, obtain continuation value.Only need step wave spectrum conversion and contravariant transducing stably to carry out the big apart from downward continuation of GEOMAGNETIC FIELD, reduce the amount of calculation of downward continuation algorithm；Recursive algorithm amount of calculation is small；Improve the big precision apart from downward continuation of the step wave number domain compensation of plane GEOMAGNETIC FIELD one.
Description
Technical field
The invention belongs to earth magnetism potential field data processing field, and in particular to a kind of plane GEOMAGNETIC FIELD downward continuation recurrence
Cosine transform method.
Background technology
In potential field data processing, magnetic anomaly downward continuation is for geologic survey, UXO detections and geomagnetic auxiliary navigation
Etc. it is significant.In sea, earth's surface and the severe area of measuring condition such as extreme low-altitude, airborne magnetic measures extremely
It is indispensable.The downward continuation of different height is carried out using Aeromagnetic data can make up inclement condition area magnetic field data not
Foot, the downward continuation of big distance can be by Aeromagnetic data and boat measurement magnetic anomaly data reduction to underwater, for the ground of submarine navigation device
Magnetic assisting navigation provides reference map.The potential field of downward continuation can protrude local magnetic anomaly, improve abnormal resolution ratio sum
According to the reliability of explanation；Downward continuation is also the core of many potential field data processing technologies, such as the calculating of normalized total gradient.
Downward continuation is a process for amplifying noise data, belongs to ill-posed problem, if prolonged using unsuitable
The method of opening up will be unable to the solution stablized.The solution that traditional downward continuation Fast Fourier Transform (FFT) method obtains is extremely unstable, can
The continuation depth leaned on be usually no more than 2~3 times of points away from.Attach most importance to apart from potential field downward continuation because its significant application value forms greatly
One of focus of magnetic data processing.Dean have studied the filtering characteristic of downward continuation from wave-number domain, utilize one group of filter coefficient
Match frequency response (the Dean W C.Frequency analysis for gravity and of downward continuation operator
magnetic interpretation, Geophysics,23(1):97-127).Clarke applications weiner equalizer designs
In the presence of downward continuation optimum filtering (the Clarke G KC.Optimum second-derivative and of interference
downward-continuation filters, Geophysics,1969,34:424-437).Ku propositions window function amendment
The filtering factor of downward continuation eliminates false anomaly (Ku C C, the Telford W formed because downward continuation amplifies High-frequency Interference
M,Lim S H.The use of linear filtering in gravity problems,Geophysics,1971,36:
1174-1203).Liang Jinwen have studied characteristic (the literary potential field downward continuations of beam brocade of four kinds of regularization downward continuation filtering factors
Regularization method, Chinese Journal of Geophysics, 1989,32 (5)： 600-608).Xu's generation Zhejiang by the potential field measured value on horizontal plane,
Initial value of the continuation face as continuation face is projected to downward vertically, with Fast Fourier Transform (FFT) upward continuation calculating observation face
Potential field value；With the potential field value on the measured value of inspection surface and the difference correction continuation face of calculated value；So iterate, until
This difference can ignore (comparison of Xu's generation Zhejiang iterative methods and FFT method potential field downward continuation effects, Chinese Journal of Geophysics, 2007,
50(1)：285-289).The first-class space iterative method in wave-number domain research Xu's generation Zhejiang in Liu east, form wave-number domain iterative method and divide
Analyse filtering characteristic, the Strict Proof convergence (wave number of the potential field downward continuations such as Liu Dongjia, Hong Tianqiu, Jia Zhihai
Domain iterative method and its convergence, Chinese Journal of Geophysics, 2009,52 (6)： 1599-1605).The introducings such as Zeng little Niu cut-off wave number
Concept, less than the downward continuation that the wave number carries out the common continuation factor, the downward of regularization continuation operator is carried out more than the wave number
Continuation (Xiaoniu Zeng, Daizhi Liu, Xihai Li, Dingxin Chen, Chao Niu.An improved
regularized downward continuation of potential field data, Journal of Applied
Geophysics,2014,106:114-118).The iterative process that Malaysian's celebrating etc. is combined using upward continuation and horizontal derivative is real
Stable downward continuation (Guoqing Ma, Cai Liu, Danian Huang, the Lili Li.A stable of existing target seeker gyro
iterative downward continuation of potential field data,Journal of Applied
Geophysics,2013,98:205–211).Gao Yuwen combinations Wavenumber domain generalized inverse algorithm and iterative method, by gradually compensating
Method realize the stable downward continuation of potential field, continuation has higher precision, and (Gao Yuwen, Luo Yao, civil and military are compensated to downward
Open up technique study and application, Chinese Journal of Geophysics, 2012,55 (8):2747-2756).
Though iterative computing method and compensating iterative method can be stablized big apart from downward continuation, higher continuation essence is obtained
Degree, but need to iterate, iterations increases the computational efficiency for having had a strong impact on downward continuation algorithm, and iteration time
Preferable downward continuation can be just approached in the case that number is more.If iterations is very few, the precision of downward continuation algorithm and sane
Property all can degradation.
The present invention proposes that a kind of the big of step wave number domain compensation of plane GEOMAGNETIC FIELD one becomes apart from downward continuation recurrence cosine
Change method.A step wave number of the downward continuation based on cosine transform is derived according to the wave-number domain iteration compensation method of downward continuation
Domain compensation operator.Using recurrence cosine transform and inverse transformation, by the step wave-number domain compensating operator and plane of vision of downward continuation
GEOMAGNETIC FIELD the GEOMAGNETIC FIELD of continuation plane is directly calculated, without iterative process.Mended compared to wave-number domain iteration
The downward continuation fourier transform method repaid, the downward continuation cosine transform method of a step wave number domain compensation need not iterate；With
Fourier transformation is compared, and cosine transform has higher energy compression performance, can reduce Gibbs boundary effects.Therefore, with
Downward continuation fourier transform method of wave-number domain iterative compensation etc. is compared, a kind of step wave number domain compensation of GEOMAGNETIC FIELD one to
Lower continuation recurrence cosine transform method, it can not only realize greatly apart from downward continuation, and the operation time of continuation algorithm and calculating
Error is reduced, and simplifies calculation process.Recurrence cosine transform also accommodates the data sequence of random length with inverse transformation
Row, and Fast Fourier Transform (FFT) usually requires that the length of data sequence is the integral number power of radix with inverse transformation, zero padding can increase
The amount of calculation of algorithm.
The content of the invention
It is an object of the invention to provide one kind to improve arithmetic accuracy, can be applied to the one of the data sequence of random length
Kind plane GEOMAGNETIC FIELD downward continuation recurrence cosine transform method.
The purpose of the present invention is realized by following technical solution：
A kind of plane GEOMAGNETIC FIELD downward continuation recurrence cosine transform method, first carry out plane of vision magnetic anomaly field measurement
The pretreatment of data, reject the measurement data of gross error；Using Interpolation Property of Radial Basis Function method by plane of vision irregularly from
The GEOMAGNETIC FIELD measurement data for dissipating distribution carries out gridding, obtains the magnetic anomaly measurement data of gridding.Recycle recurrence
The magnetic anomaly measurement data of gridding is converted to cosine transform wave spectrum by cosine transform；Utilize a step wave number of downward continuation
Domain compensation operator, obtain the magnetic anomaly wave spectrum in continuation plane；Recycle recurrence cosine inverse transformation that continuation plane is calculated
On GEOMAGNETIC FIELD.Described in comprising the following steps that, wherein for parallel computation, step 5,6 and 7 can be carried out.
Step 1, gridding parameter is chosen, plane GEOMAGNETIC FIELD measurement data is pre-processed, rejected beyond reasonable
The thick measurement error of scope.
Step 2, using Interpolation Property of Radial Basis Function method, the discrete GEOMAGNETIC FIELD measurement data on plane of vision is subjected to net
Format, obtain the gridding magnetic anomaly measurement data Δ T (n on plane of visionx,ny,z0), nx=1,2, L, Mx, ny=1,2,
L,My.Plane of vision equation is z=z0, the equation of downward continuation plane is z=z0+ Δ z, z-axis is vertically downward.
Step 3, successively the GEOMAGNETIC FIELD on plane of vision is entered every trade to row to one-dimensional recurrence cosine transform, this
The calculating order of both direction can exchange, and obtain the two-dimensional cosine transform wave spectrum Δ T of GEOMAGNETIC FIELDC(ux,uy,z0)。
One step 4, setting a small amount of ε, for characterizing the precision of downward continuation algorithm, choose damping factor α value, 0 ＜ α
＜ 1, it is 0.1 generally to take α, and the iterations n in iterative compensation method is sought according to the inequality shown in formula (1)I。
In formula
Step 5, by formula (2), by the GEOMAGNETIC FIELD wave spectrum Δ T of plane of visionC(ux,uy,z0) calculate continuation plane ground
The wave spectrum Δ T of magnetic anomaly fieldC(ux,uy,z0+Δz)。
In formula,For a step wave-number domain compensating factor of plane GEOMAGNETIC FIELD downward continuation, it is expressed
Formula is
One step wave-number domain compensating factor is equivalent to the n in wave-number domain iteration compensation methodISecondary iterative compensation calculates.
Step 6, successively to Δ TC(ux,uy,z0+ Δ z) enter every trade to row to one-dimensional recurrence anti-cosine transform, this two
The calculating order in individual direction can exchange, and obtain continuation value Δ T (x, y, the z of the GEOMAGNETIC FIELD in continuation plane0+Δz)。
The beneficial effects of the present invention are：
It is proposed a kind of step wave number domain compensation of plane GEOMAGNETIC FIELD one it is big apart from downward continuation recurrence cosine transform method
Without iteration, the advantages that amount of calculation is small, continuation distance is big and continuation precision is higher, it is big apart from downward continuation to solve compensating iterative
The problem of method needs to iterate；The present invention only needs step wave spectrum conversion and contravariant transducing stably to carry out GEOMAGNETIC FIELD
Greatly apart from downward continuation, the amount of calculation for greatly apart from the algorithm steps of downward continuation, reducing downward continuation algorithm is simplified；Position
Field wave spectral transformation and inverse transformation are small using cosine transform and the recursive algorithm of cosine inverse transformation, amount of calculation；Simultaneously compared to quick
Fourier transformation, reduce Gibbs boundary effects, improve the big distance of the step wave number domain compensation of plane GEOMAGNETIC FIELD one downwards
The precision of continuation.Recurrence cosine transform and inverse transformation are adapted to the data sequence of random length, and Fast Fourier Transform (FFT) with
Inverse transformation usually requires that the length of data sequence is the integral number power of radix, and zero padding can increase the amount of calculation of algorithm.
Brief description of the drawings
Fig. 1 is the block diagram of plane GEOMAGNETIC FIELD downward continuation recurrence cosine transform method
Fig. 2 is the flow chart of plane GEOMAGNETIC FIELD downward continuation recurrence cosine transform method
Fig. 3 is the big spheroid mould apart from downward continuation of the step wave number domain compensation based on Fourier transformation and cosine transform
Type experimental result
Fig. 4 be the step wave number domain compensation based on Fourier transformation and cosine transform big spheroid apart from downward continuation and
Cuboid mixed model experimental result
Embodiment
The embodiment of the present invention is described further below in conjunction with the accompanying drawings：
The pretreatment of step 1, plane GEOMAGNETIC FIELD measurement data.The reasonable model of setting measurement error and magnetic anomaly value
Enclose, the GEOMAGNETIC FIELD measurement data beyond this zone of reasonableness is considered as gross error, and the measurement data of gross error is given
Reject.
The gridding processing of step 2, plane GEOMAGNETIC FIELD measurement data.Gridding parameter is chosen, utilizes radial direction base letter
Number interpolation method, discrete GEOMAGNETIC FIELD measurement data is subjected to gridding, obtains the magnetic anomaly of the gridding on plane of vision
Field measurement data.
Interpolation Property of Radial Basis Function method is a kind of high interpolation method of accuracy, determines known data point in using basic function
Insert the optimal weight of grid node；Suitable for the interpolation calculation of a large amount of point datas, there is higher precision of prediction, can be preferably
Reflect the change of data.The RBF planar interpolation expression formula of interpolation point is
In formula, (x, y) be interpolation point coordinate, f0(x, y)=c1xx+c1yy+c0, (xj,yj) be sampled point coordinate, λj
For the weights corresponding to each sampled point,For RBF, | | | | it is euclideam norm.Choose thin plate spline
Function makes d as RBF·j=| | x-xj,y-yj| |, then have
The given equation constraints of sampled point is
If f (x, y) has secondary continuous derivative, then energy function is expressed as：
Energy function is minimized, the orthogonality condition for obtaining required satisfaction is
In order to solve the coefficient c in formula (1)0、c1x、c1yAnd weights λj, simultaneous formula (3) and formula (5) are obtained on coefficient
c0、c1x、c1yAnd weights λjSystem of linear equations be
The system of linear equations that coefficient matrix is symmetric positive definite matrix is solved, obtains coefficient c0、c1x、c1yAnd weights λj(j=
1,2 ..., p), substituted into formula (1) and obtain the gridding magnetic anomaly field data on inspection surface, nx=1,2 ..., Mx, ny=
1,2,…,My, the equation of plane of vision is z=z0, z-axis is vertically downward.
Step 3, using recurrence cosine transform by Δ T (nx,ny,z0) Calculation Plane GEOMAGNETIC FIELD cosine transform wave spectrum
ΔTC(ux,uy,z0).Definition
1) step presses formula (8), by Δ T (ux,uy,z0) calculate Δ T1 C(ux,ny,z0),
In formula, ux=nx, θux=(ux-1)π/Mx,WithDetermined by recurrence formula (9).
In formula, mx=1,2, L, Kx, Kx=[(Mx+ 1)/2], symbol [] represents to round.
2) step presses formula (11), by Δ T1 C(ux,ny,z0) calculate Δ TC(ux,uy,z0),
In formula, uy=ny, θuy=(uy-1)π/My,WithDetermined by recurrence formula (12).
In formula, my=1,2, L, Ky, Ky=[(My+ 1)/2], symbol [] represents to round.
Step 4, an a small amount of ε are set according to the requirement of downward continuation precision, choose damping factor α value, 0 ＜ α ＜ 1,
It is 0.1 generally to take α.Downward continuation precision is generally unsatisfactory for requiring when iterations is too small or excessive, and solution meets formula (14)
The n of shown inequalityI。
In formula
1) k=0 is walked, according to practical experience, determines nIInitial section
2) step selects sectionIntermediate pointCalculateIf meetThen stop search；
Otherwise,
3) step selects sectionIntermediate pointCalculateIf meetThen stop search；
Otherwise,
4) 2) k=k+1, a turn step continue search for step.
Step 5, by formula (15), by the cosine transform wave spectrum Δ T of inspection surface GEOMAGNETIC FIELDC(ux,uy,z0) and plane earth
One step wave-number domain compensating factor of magnetic anomaly field downward continuation directly calculates the wave spectrum Δ T of continuation plane GEOMAGNETIC FIELD fieldC
(ux,uy,z0+Δz)。
In formula,For a step wave-number domain compensating factor of plane GEOMAGNETIC FIELD downward continuation, it is expressed
Formula is
Step 6, to Δ TC(ux,uy,z0- Δ z) enter every trade to row to one-dimensional recurrence anticosine become, obtain continuation and put down
The continuation value Δ T (n of GEOMAGNETIC FIELD on facex,ny,z0+Δz)。
1) step presses formula (17), by Δ TC(ux,uy,z0+ Δ z) calculates Δ T1 C(nx,uy,z0+ Δ z),
In formula,WithDetermined by recurrence formula (18).
2) step presses formula (19), byCalculate Δ T (nx,ny,z0+ Δ z),
In formula,WithDetermined by recurrence formula (20).
Using recurrence cosine transform and inverse transformation from the big distance of step wave-number domain of GEOMAGNETIC FIELD one of plane of vision to downward
Open up, obtain the GEOMAGNETIC FIELD of continuation plane.The step wave-number domain of plane GEOMAGNETIC FIELD one is big apart from downward continuation plane earth magnetic anomaly
The block diagram of the cosine transform method of normal field is as shown in Figure 1.The step wave-number domain of plane GEOMAGNETIC FIELD one is big apart from downward continuation plane
The algorithm flow chart of the cosine transform method of GEOMAGNETIC FIELD is as shown in Figure 2.
To characterize the big continuation error apart from downward continuation plane GEOMAGNETIC FIELD of the step wave-number domain of plane GEOMAGNETIC FIELD one,
Define root-mean-square error and relative error both error criterions.
The root-mean-square error RMSE of downward continuation is carried out using Fourier transformationFFor
The root-mean-square error RMSE of downward continuation is carried out using cosine transformCFor
The relative error RE of downward continuation is carried out using Fourier transformationFFor
The relative error RE of downward continuation is carried out using cosine transformCFor
In formula, Δ T (nx,ny,z0+ Δ z) be continuation plane GEOMAGNETIC FIELD model theory value, Δ TF(nx,ny,z0+Δ
Z) it is the continuation value of the continuation plane GEOMAGNETIC FIELD obtained using Fourier transformation and inverse transformation.
Experiment one：Plane of vision is z=z0=0, the x and y coordinates scope of observation area is all [- 10000m, 10000m],
Grid number is 256 × 256, and the depth of downward continuation is Δ z=1568.6m (20 times of points away from).The position of six spheroid magnetic sources point
Not Wei [0m, 0m, 4000m], [2000m, 800m, 3400m], [1800m, -1600m, 3500m], [- 2000m, -1800m,
3600m], [- 1000m, 2000m, 3700m] and [1500m, 0m, 4300m], the intensity of magnetization is respectively 1.6 × 105A/m、1.4
×105A/m、1.3×105A/m、1.5×105A/m、1.2×105A/m and 1.8 × 105A/m, the radius of a ball be respectively 10m, 12m,
8m, 15m, 20m and 18m, magnetic dip angle distinguish 45 °, 60 °, 75 °, 45 °, 60 ° and 25 °, magnetic declination is respectively 45 °, 35 °, 15 °,
35 °, 30 °, 45 ° and 45 °.The magnetic dip angle and magnetic declination in earth's magnetic field are respectively 60 ° and 10 °.Magnetic survey noise on inspection surface is zero
Average, the high band Gaussian noise of 20nT standard deviations.
The continuation root-mean-square error of the downward continuation of sphere model data with iterations change curve as shown in figure 3,
Wherein black line is the root-mean-square error curve obtained using Fourier transformation, and red line is the root mean square obtained using cosine transform
Error curve.The continuation relative error of downward continuation is with the change curve of iterations as shown in figure 4, wherein black line is use
The relative error curve that Fourier transformation obtains, red line are the relative error curve obtained using cosine transform.As a result base is shown
It is higher than the continuation precision based on Fourier transformation and inverse transformation in the continuation precision of cosine transform and inverse transformation, single-step iteration is mended
The precision repaid is higher than the precision without iterative compensation；The error of single-step iteration backoff algorithm starts to subtract with the increase of iterations
It is small, then no longer reduce with the increase of iterations, it can also be significantly increased for relative error.
Experiment two：Two cuboid magnetic sources are added in above-mentioned emulation experiment region.The center of first cuboid is
[0m, 0m, 3000m], length and width and high respectively 80m, 40m and 40m, total magnetization intensity are 1A/m, magnetic dip angle and magnetic declination point
Wei not be 50 ° and 30 °.The center of second cuboid is [1000m, 0m, 3500m], length and width and it is high be respectively 100m,
60m and 20m, total magnetization intensity are 0.8A/m, and magnetic dip angle and magnetic declination are respectively 60 ° and 15 °.
The continuation root-mean-square error for the model data downward continuation that spheroid mixes with cuboid is bent with the change of iterations
Line is as shown in Fig. 5, and wherein black line is the root-mean-square error curve obtained using Fourier transformation, and red line is using cosine transform
Obtained root-mean-square error curve.The continuation relative error of downward continuation with iterations change curve as shown in fig. 6, its
Middle black line is the relative error curve obtained using Fourier transformation, and red line is that the relative error obtained using cosine transform is bent
Line.As a result also indicate that the continuation precision based on cosine transform and inverse transformation is higher than the continuation based on Fourier transformation and inverse transformation
Precision, the precision of single-step iteration compensation are higher than the precision without iterative compensation；The error of single-step iteration backoff algorithm starts with iteration
The increase of number and reduce, then with the increase of iterations without being obviously reduced, can also be significantly increased for relative error.
The preferred embodiments of the present invention are the foregoing is only, are not intended to limit the invention, for the skill of this area
For art personnel, the present invention can have various modifications and variations.Within the spirit and principles of the invention, that is made is any
Modification, equivalent substitution, improvement etc., should be included in the scope of the protection.
Claims (6)
1. a kind of plane GEOMAGNETIC FIELD downward continuation recurrence cosine transform method, it is characterised in that comprise the following steps
(1) pretreatment of plane GEOMAGNETIC FIELD measurement data, gridding parameter is chosen, to plane GEOMAGNETIC FIELD measurement data
Pre-processed, reject the thick measurement error beyond zone of reasonableness；
(2) the gridding processing of plane GEOMAGNETIC FIELD measurement data, using Interpolation Property of Radial Basis Function method, by plane of vision
Discrete GEOMAGNETIC FIELD measurement data carries out gridding, obtains the gridding magnetic anomaly measurement data Δ T on plane of vision
(nx,ny,z0), nx=1,2, L, Mx, ny=1,2, L, My；
(3) using recurrence cosine transform by Δ T (nx,ny,z0) Calculation Plane GEOMAGNETIC FIELD cosine transform wave spectrum Δ TC(ux,
uy,z0)；
(4) an a small amount of ε are set according to the requirement of downward continuation precision, chooses damping factor α value, 0 ＜ α ＜ 1, generally take α
For 0.1；
(5) by the cosine transform wave spectrum Δ T of inspection surface GEOMAGNETIC FIELDC(ux,uy,z0) and plane GEOMAGNETIC FIELD downward continuation
One step wave-number domain compensating factor directly calculates the wave spectrum Δ T of continuation plane GEOMAGNETIC FIELD fieldC(ux,uy,z0+Δz)；
(6) to Δ TC(ux,uy,z0- Δ z) enter every trade to row to one-dimensional recurrence anticosine become, obtain the ground in continuation plane
The continuation value Δ T (n of magnetic anomaly fieldx,ny,z0+Δz)。
A kind of 2. plane GEOMAGNETIC FIELD downward continuation recurrence cosine transform method according to claim 1, it is characterised in that
Described step (2) is specially：
The RBF planar interpolation expression formula of interpolation point is
In formula, (x, y) be interpolation point coordinate, f0(x, y)=c1xx+c1yy+c0, (xj,yj) be sampled point coordinate, λjTo be right
Should in the weights of each sampled point,For RBF, | | | | it is euclideam norm；
Thin plate spline function is chosen as RBF, makes dj=| | x-xj,y-yj| |, then have
The given equation constraints of sampled point is
If f (x, y) has secondary continuous derivative, then energy function is expressed as：
Energy function is minimized, the orthogonality condition for obtaining required satisfaction is
In order to solve the coefficient c in formula (1)0、c1x、c1yAnd weights λj, simultaneous formula (3) and formula (5) are obtained on coefficient c0、c1x、
c1yAnd weights λjSystem of linear equations be
The system of linear equations that coefficient matrix is symmetric positive definite matrix is solved, obtains coefficient c0、c1x、c1yAnd weights λj(j=1,
2 ..., p), substituted into formula (1) and obtain the gridding magnetic anomaly field data on inspection surface, nx=1,2 ..., Mx, ny=1,
2,…,My, the equation of plane of vision is z=z0, z-axis is vertically downward.
A kind of 3. plane GEOMAGNETIC FIELD downward continuation recurrence cosine transform method according to claim 1, it is characterised in that
Described step (3) is specially：
Definition
(3.1) formula (8) is pressed, by Δ T (ux,uy,z0) calculate
In formula, ux=nx, θux=(ux-1)π/Mx,WithDetermined by recurrence formula (9)；
In formula, mx=1,2, L, Kx, Kx=[(Mx+ 1)/2], symbol [] represents to round；
2) step presses formula (11), by Δ T1 C(ux,ny,z0) calculate Δ TC(ux,uy,z0),
In formula, uy=ny, θuy=(uy-1)π/My,WithDetermined by recurrence formula (12)；
In formula, my=1,2, L, Ky, Ky=[(My+ 1)/2], symbol [] represents to round；
A kind of 4. plane GEOMAGNETIC FIELD downward continuation recurrence cosine transform method according to claim 1, it is characterised in that
Described step (4) is specially：
Downward continuation precision is generally unsatisfactory for requiring when iterations is too small or excessive, and solution meets inequality shown in formula (14)
NI
In formula
(4.1) k=0, according to practical experience, n is determinedIInitial section
(4.2) section is selectedIntermediate pointCalculateIf meetThen stop search；It is no
Then,
(4.3) section is selectedIntermediate pointCalculateIf meetThen stop search；It is no
Then,
(4.4) k=k+1, turn step and 2) continue search for.
A kind of 5. plane GEOMAGNETIC FIELD downward continuation recurrence cosine transform method according to claim 1, it is characterised in that
Described step (5) is specially：
In formula,For a step wave-number domain compensating factor of plane GEOMAGNETIC FIELD downward continuation, its expression formula is
A kind of 6. plane GEOMAGNETIC FIELD downward continuation recurrence cosine transform method according to claim 1, it is characterised in that
Described step (6) is specially：
(6.1) formula (17) is pressed, by Δ TC(ux,uy,z0+ Δ z) is calculated
In formula,WithDetermined by recurrence formula (18)；
(6.2) formula (19) is pressed, by Δ T1 C(nx,uy,z0+ Δ z) calculates Δ T (nx,ny,z0+ Δ z),
In formula,WithDetermined by recurrence formula (20)；
The root-mean-square error RMSE of downward continuation is carried out using Fourier transformationFFor
The root-mean-square error RMSE of downward continuation is carried out using cosine transformCFor
The relative error RE of downward continuation is carried out using Fourier transformationFFor
The relative error RE of downward continuation is carried out using cosine transformCFor
In formula, Δ T (nx,ny,z0+ Δ z) be continuation plane GEOMAGNETIC FIELD model theory value, Δ TF(nx,ny,z0+ Δ z) is
The continuation value of the continuation plane GEOMAGNETIC FIELD obtained using Fourier transformation and inverse transformation.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710795560.4A CN107632964B (en) | 2017-09-06 | 2017-09-06 | Downward continuation recursive cosine transform method for plane geomagnetic abnormal field |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710795560.4A CN107632964B (en) | 2017-09-06 | 2017-09-06 | Downward continuation recursive cosine transform method for plane geomagnetic abnormal field |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN107632964A true CN107632964A (en) | 2018-01-26 |
| CN107632964B CN107632964B (en) | 2021-06-11 |
Family
ID=61100103
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201710795560.4A Active CN107632964B (en) | 2017-09-06 | 2017-09-06 | Downward continuation recursive cosine transform method for plane geomagnetic abnormal field |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN107632964B (en) |
Cited By (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN109190082A (en) * | 2018-08-15 | 2019-01-11 | 中国人民解放军61540部队 | Inhibit the cosine function flared end method of boundary effect when airborne gravity data downward continuation |
| CN110389391A (en) * | 2019-08-01 | 2019-10-29 | 自然资源部第二海洋研究所 | A kind of heavy magnetic potential field desorptiion continuation method based on spatial domain |
| CN110531430A (en) * | 2019-08-29 | 2019-12-03 | 中国石油天然气集团公司 | Processing method, device and the electronic equipment of submarine pipeline magnetic survey data |
| CN110543611A (en) * | 2019-08-15 | 2019-12-06 | 桂林理工大学 | Low latitude magnetic anomaly data magnetization pole calculation method and device |
| CN111337992A (en) * | 2020-03-23 | 2020-06-26 | 兰州大学 | Method for obtaining depth of field source based on downward continuation of bit field data |
| CN114048636A (en) * | 2022-01-11 | 2022-02-15 | 中国测绘科学研究院 | Gravity anomaly calculation method and device based on wavelet transformation |
| CN114624781A (en) * | 2022-03-17 | 2022-06-14 | 宿州学院 | Method and device for extracting ocean magnetic abnormal strip under magma activity interference |
| CN116774301A (en) * | 2023-06-26 | 2023-09-19 | 中国人民解放军火箭军工程大学 | Method, system, electronic equipment and medium for downward continuation of heavy magnetic bit field regularization |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20080052562A1 (en) * | 2006-08-22 | 2008-02-28 | Denso Corporation | Fault detection unit for rotation angle detecting device |
| CN103900564A (en) * | 2014-03-03 | 2014-07-02 | 哈尔滨工程大学 | Submergence assisted geomagnetic anomaly inversion velocity measurement/underwater continuous positioning method |
| CN104215259A (en) * | 2014-08-22 | 2014-12-17 | 哈尔滨工程大学 | Inertial navigation error correction method based on geomagnetism modulus gradient and particle filter |
-
2017
- 2017-09-06 CN CN201710795560.4A patent/CN107632964B/en active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20080052562A1 (en) * | 2006-08-22 | 2008-02-28 | Denso Corporation | Fault detection unit for rotation angle detecting device |
| CN103900564A (en) * | 2014-03-03 | 2014-07-02 | 哈尔滨工程大学 | Submergence assisted geomagnetic anomaly inversion velocity measurement/underwater continuous positioning method |
| CN104215259A (en) * | 2014-08-22 | 2014-12-17 | 哈尔滨工程大学 | Inertial navigation error correction method based on geomagnetism modulus gradient and particle filter |
Non-Patent Citations (8)
| Title |
|---|
| HUI ZHANG ET AL: "Analysis on Convergence of Iteration Method for Potential Fields Downward Continuation and Research on Robust Downward Continuation Method", 《CHINESE JOURNAL OF GEOPHYSICS》 * |
| SHI‐ZHE XU ET AL: "Downward Continuation of Aeromagnetic Anomaly from Flying Altitude to Terrain", 《CHINESE JOURNAL OF GEOPHYSICS》 * |
| 张凤旭 等: "利用余弦变换计算重力异常的向上延拓", 《地球物理学进展》 * |
| 张凤琴: "基于离散余弦变换的位场谱方法及应用", 《中国博士学位论文全文数据库基础科学辑》 * |
| 曾小牛 等: "位场向下延拓三种迭代方法之比较", 《地球物理学进展》 * |
| 沈宏君: "任意长离散余弦变换的快速递归算法", 《电子与信息学报》 * |
| 王宝龙 等，: "正则化薄板样条函数拟合地层界面", 《第十五届全国数学地质与地学信息学术研讨会》 * |
| 马涛: "地磁异常数据向下延拓算法研究", 《中国优秀硕士学位论文全文数据库信息科技辑》 * |
Cited By (13)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN109190082B (en) * | 2018-08-15 | 2022-12-30 | 中国人民解放军61540部队 | Method for downwardly extending aviation gravity data |
| CN109190082A (en) * | 2018-08-15 | 2019-01-11 | 中国人民解放军61540部队 | Inhibit the cosine function flared end method of boundary effect when airborne gravity data downward continuation |
| CN110389391B (en) * | 2019-08-01 | 2020-12-15 | 自然资源部第二海洋研究所 | Heavy magnetic bit field analytic extension method based on spatial domain |
| CN110389391A (en) * | 2019-08-01 | 2019-10-29 | 自然资源部第二海洋研究所 | A kind of heavy magnetic potential field desorptiion continuation method based on spatial domain |
| CN110543611B (en) * | 2019-08-15 | 2022-11-25 | 桂林理工大学 | Low latitude magnetic abnormal data magnetization pole calculation method and device |
| CN110543611A (en) * | 2019-08-15 | 2019-12-06 | 桂林理工大学 | Low latitude magnetic anomaly data magnetization pole calculation method and device |
| CN110531430A (en) * | 2019-08-29 | 2019-12-03 | 中国石油天然气集团公司 | Processing method, device and the electronic equipment of submarine pipeline magnetic survey data |
| CN111337992A (en) * | 2020-03-23 | 2020-06-26 | 兰州大学 | Method for obtaining depth of field source based on downward continuation of bit field data |
| CN111337992B (en) * | 2020-03-23 | 2021-04-06 | 兰州大学 | Method for obtaining depth of field source based on downward continuation of bit field data |
| CN114048636A (en) * | 2022-01-11 | 2022-02-15 | 中国测绘科学研究院 | Gravity anomaly calculation method and device based on wavelet transformation |
| CN114624781A (en) * | 2022-03-17 | 2022-06-14 | 宿州学院 | Method and device for extracting ocean magnetic abnormal strip under magma activity interference |
| CN116774301A (en) * | 2023-06-26 | 2023-09-19 | 中国人民解放军火箭军工程大学 | Method, system, electronic equipment and medium for downward continuation of heavy magnetic bit field regularization |
| CN116774301B (en) * | 2023-06-26 | 2024-05-14 | 中国人民解放军火箭军工程大学 | Method, system, electronic equipment and medium for downward continuation of heavy magnetic bit field regularization |
Also Published As
| Publication number | Publication date |
|---|---|
| CN107632964B (en) | 2021-06-11 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN107632964A (en) | A kind of plane GEOMAGNETIC FIELD downward continuation recurrence cosine transform method | |
| CN107230197B (en) | Tropical cyclone objective strength determination method based on satellite cloud image and RVM | |
| CN107291659A (en) | The recurrence cosine transform method of the step upward continuation plane modulus gradient of plane GEOMAGNETIC FIELD one | |
| CN104268880A (en) | Depth information obtaining method based on combination of features and region matching | |
| CN107607195B (en) | A kind of beam quality measurement method obtained in real time based on complex amplitude | |
| CN104867179A (en) | Whole spectral range optical imager remote sensing image simulation method | |
| CN103364833A (en) | High-precision dip estimation method | |
| Liu et al. | NLOS signal detection and correction for smartphone using convolutional neural network and variational mode decomposition in urban environment | |
| Zhang et al. | A study on coastline extraction and its trend based on remote sensing image data mining | |
| CN109738852A (en) | The distributed source two-dimensional space Power estimation method rebuild based on low-rank matrix | |
| CN106323271A (en) | Spacecraft relative attitude measurement vector selection method based on feature singular values | |
| Wang et al. | Interference image registration combined by enhanced scale-invariant feature transform characteristics and correlation coefficient | |
| CN105225232A (en) | A kind of colour of view-based access control model attention mechanism warship cooperative target detection method | |
| Liu et al. | Laser spot center location algorithm based on sub-pixel interpolation | |
| CN105096328A (en) | Camera self-calibration method based on IOS-PSO (Internetwork Operating System-Particle Swarm Optimization) | |
| Wu et al. | Noniterative calibration of a camera lens with radial distortion | |
| Xu | [Retracted] Application of Remote Sensing Image Data Scene Generation Method in Smart City | |
| Vitkin et al. | Study of in situ calibration performance of co-located multi-sensor hot-film and sonic anemometers using a ‘virtual probe’algorithm | |
| Avitzour | Novel scene calibration procedure for video surveillance systems | |
| CN114488217A (en) | High-orbit satellite CEI signal frequency estimation method based on deep learning | |
| Rosenbaum | Enhanced propagation modeling of directional aviation noise: A hybrid parabolic equation-fast field program method | |
| CN105628342B (en) | A kind of optical flow visualization technical optics are laid out method of value solving | |
| CN111044489A (en) | Method for obtaining atmosphere refractive index height distribution profile based on multi-wavelength measurement | |
| CN114239379B (en) | Power transmission line geological disaster analysis method and system based on deformation detection | |
| Arnold et al. | Vision-Assisted Digital Twin Creation for mmWave Beam Management |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |