CN106197405A - Inertia earth magnetism matching locating method under the influence of geomagnetic diurnal change - Google Patents

Abstract

The invention provides the inertia earth magnetism matching locating method under the influence of a kind of geomagnetic diurnal change, read the position measurements a of the point to be matched in N number of moment from inertial navigation system_iAnd b_i, obtain geomagnetic field intensity information I from gaussmeter_i；The position of the point N number of to be matched according to inertial navigation system instruction, reads the reference value I (a of corresponding geomagnetic field intensity from the geomagnetic database prestored respectively_i,b_i), reference value I of the gradient of geomagnetic field intensity_x,iAnd I_y,i；Introduce and initialize longitude and latitude error, course error and geomagnetic diurnal change error；Calculate the increment of longitude and latitude error, the increment of course error and the increment δ M of geomagnetic diurnal change error；Update longitude and latitude error, course error and geomagnetic diurnal change error M；Judge whether to meet and terminate iterated conditional, calculate parameter K and δ K according to the M after updating；Obtain longitude and latitude error, course error and geomagnetic diurnal change error according to iterative computation, acquired results substitution coupling track is i.e. obtained with the relation equation of reference locus and mates track.

Description

Inertia earth magnetism matching locating method under the influence of geomagnetic diurnal change

Technical field

The present invention relates to the inertia earth magnetism matching locating method under the influence of a kind of geomagnetic diurnal change, belong to inertia geomagnetic matching fixed Method for position technical field.

Background technology

Airmanship along with modern science and technology development make rapid progress, mainly have inertial navigation, satellite navigation, The airmanships such as shape coupling navigation, celestial navigation, earth-magnetism navigation.In airmanship miscellaneous, earth-magnetism navigation is high with it Degree autonomy, disguise and become the one of current navigation field without advantages such as accumulated errors and study greatly hot topic.Inertia geomagnetic matching Referring to, carrier is provided with inertial navigation system and Magnetic Sensor, carrier is in flight course, and inertial navigation system exports the motion rail of carrier Mark information, is referred to as " reference locus " by this output trajectory, owing to inertia device exists drift, therefore refers to track and carrier movement There is certain error in the real trace experienced；While inertial navigation system output carrier track, measure geomagnetic sensor To the geomagnetic field information of the experienced position of carrier, magnetic field of the earth information, earth's magnetic field that then basis is pre-stored in computer are real-time Metrical information and carrier reference locus, use matching process, obtain the coupling track of carrier.Earth's magnetic field comprise stabilizing magnetic field and Variation magnetic field two parts.The former includes main field and anomalous field, is the main body in earth's magnetic field, is also geomagnetic matching letter to be used Breath；The latter originates from the space current system outside solid earth, is equivalent to interference during magnetic survey.After geomagnetic diurnal change belongs to Person, changes including solar quiet day change and lunar day, and its mean range is about a few nanotesla to tens nanoteslas.At present, earth magnetism Day becomes and does not form the most unified mathematical model, and the geomagnetic diurnal change curve of different location has bigger difference, same place Geomagnetic diurnal change curve Various Seasonal is the most different in also different, daytime and evening, is therefore difficult in inertia geomagnetic matching position fixing process Effectively compensate or reject.So, the precision of inertia geomagnetic matching location and reliability can be produced tight by geomagnetic diurnal change The impact of weight.

For these reasons, the most urgent to the needs of the inertia earth magnetism matching locating method under the influence of geomagnetic diurnal change. Yet with the uncertainty become day, complexity, day is become and is extremely difficult to process, to this respect in existing earth-magnetism navigation method Process considerably less.The document having in current earth-magnetism navigation research utilizes complicated approach matching geomagnetic diurnal change curve, has Although document emphasizes that geomagnetic diurnal change affects, but could not enough provide the method processing geomagnetic diurnal change, directly ignore geomagnetic diurnal change Impact.By literature search, Authorization Notice No. be CN103115624B, authorized announcement date be in December, 2014 Chinese invention of 10 days Patent, discloses " a kind of geomagnetic diurnal change modification method based on geomagnetic matching ", and day change is carried out by this patent based on FMI method Matching, same day and calculating based on the geomagnetic data of the most totally three days and the method needs to navigate, be difficult to actual should With.Xie Shimin etc. indicate that GEOMAGNETIC CHANGE field (namely geomagnetic diurnal change) is right in article " geomagnetic matching navigation Key Technology Analysis " The navigation accuracy bigger error of generation, but the method that can not propose to correct geomagnetic diurnal change.Li Yuze etc. at article " based on ICCP The geomagnetic matching localization method of algorithm " (see " scientific algorithm and information processing ") does not the most mention geomagnetic diurnal change, have ignored ground The impact that magnetic day becomes, it is impossible to ensure coupling positioning precision.

Summary of the invention

The invention aims to solve the problem that above-mentioned prior art exists, and then a kind of geomagnetic diurnal change impact is provided Under inertia earth magnetism matching locating method.

It is an object of the invention to be achieved through the following technical solutions:

A kind of inertia earth magnetism matching locating method under the influence of geomagnetic diurnal change, step is as follows:

Step one, read current time and the position measurements a of the point to be matched in front N-1 moment from inertial navigation system_i And b_i, wherein a_iRepresent longitude, b_iRepresenting latitude, subscript i represents the most in the same time, i=1 ... N, N are integer and N ＞ 2, and i is 1 table Show current time, obtain current time and the measured value I of the geomagnetic field intensity in front N-1 moment by geomagnetic sensor_i；

Step 2, the position of the point N number of to be matched indicated according to inertial navigation system, divide from the geomagnetic database prestored Do not read the reference value I (a of corresponding geomagnetic field intensity_i,b_i), reference value I of the gradient of geomagnetic field intensity_x,iAnd I_y,i, wherein I_x,iRepresent the geomagnetic field intensity gradient in longitudinal value on i-th position, I_y,iRepresent that geomagnetic field intensity is at latitude The gradient in direction value on i-th position；

Step 3, introduce and initialize longitude and latitude error, course error and geomagnetic diurnal change error:

M=[Δ x Δ y α δ]^T=[0 00 0]^T

Step 4, calculate iterative parameter G, F and H according to formula (1), formula (2) and formula (3):

G=g (M) (1)

F=f (M) (2)

H=F^-1 (3)

Wherein:

G (M)=[g₁(M) g₂(M) g₃(M) g₄(M)]^T；

$f\left(M\right)=\left[\begin{array}{cccc}{f}_{11}\left(M\right)& f_{12}\left(M\right)& f_{13}\left(M\right)& f_{14}\left(M\right)\\ f_{12}\left(M\right)& f_{22}\left(M\right)& f_{23}\left(M\right)& f_{24}\left(M\right)\\ f_{13}\left(M\right)& f_{23}\left(M\right)& f_{33}\left(M\right)& f_{34}\left(M\right)\\ f_{14}\left(M\right)& f_{24}\left(M\right)& f_{34}\left(M\right)& f_{44}\left(M\right)\end{array}\right];$

$\begin{array}{c}{g}_1\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{x,i}\[I_{x,i}(a_i^{\′}\cos \mathrm{\α}+b_i^{\′}\sin \mathrm{\α}-a_i^{\′}+\mathrm{\Δ}x)\\ +I_{y,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α}-b_i^{\′}+\mathrm{\Δ}y)+I_{t,i}+\mathrm{\δ}\]\end{array};$

$\begin{array}{c}{g}_2\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{y,i}\[I_{x,i}(a_i^{\′}\cos \mathrm{\α}+b_i^{\′}\sin \mathrm{\α}-a_i^{\′}+\mathrm{\Δ}x)\\ +I_{y,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α}-b_i^{\′}+\mathrm{\Δ}y)+I_{t,i}+\mathrm{\δ}\]\end{array};$

$\begin{array}{c}{g}_3\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}\[I_{x,i}(a_i^{\′}\cos \mathrm{\α}+b_i^{\′}\sin \mathrm{\α}-a_i^{\′}+\mathrm{\Δ}x)+I_{y,i}(-a_i^{\′}\sin \mathrm{\α}\\ +b_i^{\′}\cos \mathrm{\α}-b_i^{\′}+\mathrm{\Δ}y)+I_{t,i}+\mathrm{\δ}\]\×\[I_{x,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α})\\ +I_{y,i}(-a_i^{\′}\cos \mathrm{\α}-b_i^{\′}\sin \mathrm{\α})\]\end{array};$

$\begin{array}{c}{g}_4\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}\[I_{x,i}(a_i^{\′}\cos \mathrm{\α}+b_i^{\′}\sin \mathrm{\α}-a_i^{\′}+\mathrm{\Δ}x)\\ +I_{y,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α}-b_i^{\′}+\mathrm{\Δ}y)+I_{t,i}+\mathrm{\δ}\]\end{array};$

$f_{11}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{x,i}^2;$

$f_{12}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{x,i}{I}_{y,i};$

$f_{13}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{x,i}\[I_{x,i}(-a_i^{\′}sin\mathrm{\α}+b_i^{\′}cos\mathrm{\α})+I_{y,i}(-a_i^{\′}cos\mathrm{\α}-b_i^{\′}sin\mathrm{\α})\];$

$f_{14}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{x,i};$

$f_{22}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{y,i}^2;$

$f_{23}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{y,i}\[I_{x,i}(-a_i^{\′}sin\mathrm{\α}+b_i^{\′}cos\mathrm{\α})+I_{y,i}(-a_i^{\′}cos\mathrm{\α}-b_i^{\′}sin\mathrm{\α})\];$

$f_{24}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{y,i};$

$\begin{array}{c}{f}_{33}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{\[I_{x,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α})+I_{y,i}(-a_i^{\′}\cos \mathrm{\α}-b_i^{\′}\sin \mathrm{\α})\]}^2\\ +\[I_{x,i}(a_i^{\′}\cos \mathrm{\α}+b_i^{\′}\sin \mathrm{\α}-a_i^{\′}+\mathrm{\Δ}x)+I_{y,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α}\\ -b_i^{\′}+\mathrm{\Δ}y)+I_{t,i}+\mathrm{\δ}\]\×\[I_{x,i}(-a_i^{\′}\cos \mathrm{\α}-b_i^{\′}\sin \mathrm{\α})+I_{y,i}(a_i^{\′}\sin \mathrm{\α}\\ -b_i^{\′}\cos \mathrm{\α})\]\end{array};$

$f_{34}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}\[I_{x,i}(-a_i^{\′}sin\mathrm{\α}+b_i^{\′}cos\mathrm{\α})+I_{y,i}(-a_i^{\′}cos\mathrm{\α}-b_i^{\′}sin\mathrm{\α})\];$

f₄₄(M)=N；

a_i'=a_i-a₁,b_i'=b_i-b₁,I_t,i=I (a_i,b_i)-I_i；

Step 5, the calculating increment of longitude and latitude error, the increment of course error and the increment δ M of geomagnetic diurnal change error:

δ M=-H × G (4)

Step 6, renewal longitude and latitude error, course error and geomagnetic diurnal change error M:

M=M+ δ M (5)

Step 7, judge whether to meet and terminate iterated conditional, if meeting, stopping iteration and also jumping to step 10, otherwise jumping to Step 8；

Terminate iterated conditional for 1., 2. in any one or two: 1. iterations reaches preset times；2. longitude and latitude is by mistake 2 norms of the increment δ M of increment, the increment of course error and the geomagnetic diurnal change error of difference are less than setting value, i.e.

||δM||₂＜ ε (6)

Wherein ε is iteration minimum error set in advance；

Step 8, calculate parameter K and δ K according to the M after updating:

K=g (M) (7)

δ K=K-G (8)

Step 9, update iteration variable H and G according to formula (9) and formula (10), then jump to step 5,

H=H+ [(δ M)-H (δ K)] (δ M)^TH/[(δM)^TH(δK)] (9)

G=K (10)

Step 10, the longitude and latitude error that obtains according to iterative computation, course error and geomagnetic diurnal change error delta x, Δ y, α and δ, Acquired results substitution coupling track is i.e. obtained with the relation equation (11) of reference locus and mates track；

$\left[\begin{array}{c}{u}_i\\ v_i\end{array}\right]=\left[\begin{array}{cc}cos\mathrm{\α}& sin\mathrm{\α}\\ -sin\mathrm{\α}& cos\mathrm{\α}\end{array}\right]\left[\begin{array}{c}{a}_i-a_1\\ b_i-b_1\end{array}\right]+\left[\begin{array}{c}{a}_1\\ b_1\end{array}\right]+\left[\begin{array}{c}\mathrm{\Δ}x\\ \mathrm{\Δ}y\end{array}\right]---\left(11\right)$

Wherein u_iIt is the position longitude of the i-th moment matching result, v_iIt it is the position latitude of the i-th moment matching result.

Geomagnetic chart in reference locus that the present invention exports using inertial navigation system, computer, geomagnetic sensor measured value as Input, calculates reference locus initial longitude and latitude error, course with the inertia earth magnetism matching locating method under the influence of geomagnetic diurnal change Error, longitude and latitude error resolving obtained, course error substitute into the relation equation of coupling track and reference locus and i.e. obtain Join result.Inertia earth magnetism matching locating method under the influence of geomagnetic diurnal change proposed by the invention has higher positioning precision, And after iteration, positioning precision has obtained further raising.

Accompanying drawing explanation

Fig. 1 is the flow chart implementing the inventive method as a example by certain experiment sport car.

Fig. 2 is longitude error curve chart.

Fig. 3 is latitude error curve chart.

Detailed description of the invention

The present invention is described in further detail below: the present embodiment enters under premised on technical solution of the present invention Row is implemented, and gives detailed embodiment, but protection scope of the present invention is not limited to following embodiment.

Inertia earth magnetism matching locating method under the influence of a kind of geomagnetic diurnal change involved by the present embodiment, including following step Rapid:

Step one, read current time and the position measurements a of the point to be matched in front N-1 moment from inertial navigation system_i And b_i, wherein a_iRepresent longitude, b_iRepresenting latitude, subscript i represents the most in the same time, i=1 ... N, N are integer and N ＞ 2, and i is 1 table Show current time, obtain current time and the measured value I of the geomagnetic field intensity in front N-1 moment by geomagnetic sensor_i；

Step 2, the position of the point N number of to be matched indicated according to inertial navigation system, divide from the geomagnetic database prestored Do not read the reference value I (a of corresponding geomagnetic field intensity_i,b_i), reference value I of the gradient of geomagnetic field intensity_x,iAnd I_y,i, wherein I_x,iRepresent the geomagnetic field intensity gradient in longitudinal value on i-th position, I_y,iRepresent that geomagnetic field intensity is at latitude The gradient in direction value on i-th position；

Step 3, introduce and initialize longitude and latitude error, course error and geomagnetic diurnal change error:

M=[Δ x Δ y α δ]^T=[0 00 0]^T

Step 4, calculate iterative parameter G, F and H according to formula (1), formula (2) and formula (3):

G=g (M) (1)

F=f (M) (2)

H=F^-1 (3)

Wherein:

G (M)=[g₁(M) g₂(M) g₃(M) g₄(M)]^T；

$f\left(M\right)=\left[\begin{array}{cccc}{f}_{11}\left(M\right)& f_{12}\left(M\right)& f_{13}\left(M\right)& f_{14}\left(M\right)\\ f_{12}\left(M\right)& f_{22}\left(M\right)& f_{23}\left(M\right)& f_{24}\left(M\right)\\ f_{13}\left(M\right)& f_{23}\left(M\right)& f_{33}\left(M\right)& f_{34}\left(M\right)\\ f_{14}\left(M\right)& f_{24}\left(M\right)& f_{34}\left(M\right)& f_{44}\left(M\right)\end{array}\right];$

$\begin{array}{c}{g}_1\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{x,i}\[I_{x,i}(a_i^{\′}\cos \mathrm{\α}+b_i^{\′}\sin \mathrm{\α}-a_i^{\′}+\mathrm{\Δ}x)\\ +I_{y,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α}-b_i^{\′}+\mathrm{\Δ}y)+I_{t,i}+\mathrm{\δ}\]\end{array};$

$\begin{array}{c}{g}_2\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{y,i}\[I_{x,i}(a_i^{\′}\cos \mathrm{\α}+b_i^{\′}\sin \mathrm{\α}-a_i^{\′}+\mathrm{\Δ}x)\\ +I_{y,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α}-b_i^{\′}+\mathrm{\Δ}y)+I_{t,i}+\mathrm{\δ}\]\end{array};$

$\begin{array}{c}{g}_3\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}\[I_{x,i}(a_i^{\′}\cos \mathrm{\α}+b_i^{\′}\sin \mathrm{\α}-a_i^{\′}+\mathrm{\Δ}x)+I_{y,i}(-a_i^{\′}\sin \mathrm{\α}\\ +b_i^{\′}\cos \mathrm{\α}-b_i^{\′}+\mathrm{\Δ}y)+I_{t,i}+\mathrm{\δ}\]\×\[I_{x,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α})\\ +I_{y,i}(-a_i^{\′}\cos \mathrm{\α}-b_i^{\′}\sin \mathrm{\α})\]\end{array};$

$\begin{array}{c}{g}_4\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}\[I_{x,i}(a_i^{\′}\cos \mathrm{\α}+b_i^{\′}\sin \mathrm{\α}-a_i^{\′}+\mathrm{\Δ}x)\\ +I_{y,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α}-b_i^{\′}+\mathrm{\Δ}y)+I_{t,i}+\mathrm{\δ}\]\end{array};$

$f_{11}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{x,i}^2;$

$f_{12}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{x,i}{I}_{y,i};$

$f_{13}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{x,i}\[I_{x,i}(-a_i^{\′}sin\mathrm{\α}+b_i^{\′}cos\mathrm{\α})+I_{y,i}(-a_i^{\′}cos\mathrm{\α}-b_i^{\′}sin\mathrm{\α})\];$

$f_{14}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{x,i};$

$f_{22}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{y,i}^2;$

$f_{23}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{y,i}\[I_{x,i}(-a_i^{\′}sin\mathrm{\α}+b_i^{\′}cos\mathrm{\α})+I_{y,i}(-a_i^{\′}cos\mathrm{\α}-b_i^{\′}sin\mathrm{\α})\];$

$f_{24}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{y,i};$

$\begin{array}{c}{f}_{33}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{\[I_{x,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α})+I_{y,i}(-a_i^{\′}\cos \mathrm{\α}-b_i^{\′}\sin \mathrm{\α})\]}^2\\ +\[I_{x,i}(a_i^{\′}\cos \mathrm{\α}+b_i^{\′}\sin \mathrm{\α}-a_i^{\′}+\mathrm{\Δ}x)+I_{y,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α}\\ -b_i^{\′}+\mathrm{\Δ}y)+I_{t,i}+\mathrm{\δ}\]\×\[I_{x,i}(-a_i^{\′}\cos \mathrm{\α}-b_i^{\′}\sin \mathrm{\α})+I_{y,i}(a_i^{\′}\sin \mathrm{\α}\\ -b_i^{\′}\cos \mathrm{\α})\]\end{array};$

$f_{34}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}\[I_{x,i}(-a_i^{\′}sin\mathrm{\α}+b_i^{\′}cos\mathrm{\α})+I_{y,i}(-a_i^{\′}cos\mathrm{\α}-b_i^{\′}sin\mathrm{\α})\];$

f₄₄(M)=N；

a_i'=a_i-a₁,b_i'=b_i-b₁,I_t,i=I (a_i,b_i)-I_i；

Step 5, the calculating increment of longitude and latitude error, the increment of course error and the increment δ M of geomagnetic diurnal change error:

δ M=-H × G (4)

Step 6, renewal longitude and latitude error, course error and geomagnetic diurnal change error M:

M=M+ δ M (5)

Step 7, judge whether to meet and terminate iterated conditional, if meeting, stopping iteration and also jumping to step 10, otherwise jumping to Step 8；

Terminate iterated conditional for 1., 2. in any one or two: 1. iterations reaches preset times；2. longitude and latitude is by mistake 2 norms of the increment δ M of increment, the increment of course error and the geomagnetic diurnal change error of difference are less than setting value, i.e.

||δM||₂＜ ε (6)

Wherein ε is iteration minimum error set in advance；

Step 8, calculate parameter K and δ K according to the M after updating:

K=g (M) (7)

δ K=K-G (8)

Step 9, update iteration variable H and G according to formula (9) and formula (10), then jump to step 5,

H=H+ [(δ M)-H (δ K)] (δ M)^TH/[(δM)^TH(δK)] (9)

G=K (10)

Step 10, the longitude and latitude error that obtains according to iterative computation, course error and geomagnetic diurnal change error delta x, Δ y, α and δ, Acquired results substitution coupling track is i.e. obtained with the relation equation (11) of reference locus and mates track；

$\left[\begin{array}{c}{u}_i\\ v_i\end{array}\right]=\left[\begin{array}{cc}cos\mathrm{\α}& sin\mathrm{\α}\\ -sin\mathrm{\α}& cos\mathrm{\α}\end{array}\right]\left[\begin{array}{c}{a}_i-a_1\\ b_i-b_1\end{array}\right]+\left[\begin{array}{c}{a}_1\\ b_1\end{array}\right]+\left[\begin{array}{c}\mathrm{\Δ}x\\ \mathrm{\Δ}y\end{array}\right]---\left(11\right)$

Wherein u_iIt is the position longitude of the i-th moment matching result, v_iIt it is the position latitude of the i-th moment matching result.

Geomagnetic field intensity described previously can select Geomagnetic Total Field, GEOMAGNETIC FIELD overall strength or earth's magnetic field the strongest Degree component of a direction under geographic coordinate system；When geomagnetic field intensity is Geomagnetic Total Field, gaussmeter uses scalar magnetic Strong meter or three axial vector gaussmeters, directly obtain Geomagnetic Total Field from described gaussmeter, as the measurement of geomagnetic field intensity Value I_i；Correspondingly, the geomagnetic field intensity prestored and the gradient of geomagnetic field intensity should be Geomagnetic Total Field and earth's magnetic field is total The gradient of intensity.When geomagnetic field intensity is GEOMAGNETIC FIELD overall strength, gaussmeter uses scalar magnetometer or three axial vectors Gaussmeter, directly obtains Geomagnetic Total Field from gaussmeter, and base area signal magnetic field model calculates GEOMAGNETIC FIELD overall strength, Measured value I as geomagnetic field intensity_i；Correspondingly, the geomagnetic field intensity prestored and the gradient of geomagnetic field intensity should be ground Magnetic anomaly field overall strength and the gradient of GEOMAGNETIC FIELD overall strength.When geomagnetic field intensity is that Geomagnetic Total Field is in geographic coordinate system During the component of lower a direction, gaussmeter uses three axial vector gaussmeters, according to three axial vector measured value and loads of gaussmeter The attitude of body, calculates Geomagnetic Total Field component of the direction under geographic coordinate system, as the measured value of geomagnetic field intensity I_i；Correspondingly, the geomagnetic field intensity prestored and the gradient of geomagnetic field intensity should be Geomagnetic Total Field in geographic coordinate system The component of the lower direction and the gradient of this component.

Inertia earth magnetism matching locating method under the influence of geomagnetic diurnal change can also include second iteration step, and described secondary is repeatedly Rapid content of riding instead of walk is: second time performs described step one to ten, with perform for the first time (i.e. an iteration) have following 2 not With: the coupling positioning result u obtained by 1. performing with first time_iAnd v_iSubstitute a read in described step one from inertial navigation system_i And b_i(as shown in formula (12)), will initial as second iteration of revised to be matched some positional value obtaining of an iteration Value；2. use the geomagnetic diurnal change error delta obtained by performing for the first time as the initial geomagnetic diurnal change error in described step 3 (such as formula (13) shown in).Owing to the initial value of second iteration is to have eliminated major part initial position through what first time iteration obtained Error, the value of major part initial heading error, again by iterative computation, its result ratio iteration for the first time is closer to true rail Mark, makes coupling positioning result precision be further enhanced.

$\left\{\begin{array}{c}{a}_i=u_i\\ b_i=v_i\end{array}\right.---\left(12\right)$

M=[0 00 δ]^T (13)

According to Fig. 1, as a example by certain experiment sport car, implement the inertia geomagnetic matching location side under the influence of geomagnetic diurnal change of the present invention The process of method is as follows:

Experiment condition: earth's magnetic field optionally magnetic anomaly field overall strength.Proton magnetometer is selected to measure Magnetic Field, matter in real time The main performance index of sub-magnetometer is as follows: resolution: 0.01nT, precision: ± 0.2nT.The main performance index of inertial navigation system As follows: the inclined unstability of gyro zero: 0.01 °/h, Gyro Random migration:Accelerometer bias unstability: 80 μ G, accelerometer random walk:

GEOMAGNETIC FIELD overall strength data interval to longitude range 126 ° to 129 ° and latitude scope 44 ° to 46 ° are stored in Airborne computer, uses forward difference method calculate the gradient information of GEOMAGNETIC FIELD overall strength and be stored in airborne computer；Choosing Taking to be matched some number is 7, i.e. N=7.

Use the inertia earth magnetism matching locating method step under the influence of geomagnetic diurnal change as follows:

Step one: read current time and the position measurements of the point to be matched in front 6 moment from inertial navigation system a_iAnd b_i, as shown in table 1；Measured value according to proton magnetometer and magnetic field of the earth model, obtain current time and when first 6 Magnetic field intensity metrical information I carved_iAs shown in table 2:

Table 1 inertial navigation system measures position coordinates

i	Longitude ai(°)	Latitude bi(°)
			1	127.1558	44.3980
2	127.1988	44.4154
			3	127.2619	44.4456
4	127.3237	44.4746
			5	127.3836	44.5033
6	127.4410	44.5324
			7	127.4946	44.5630

Table 2 magnetic field intensity metrical information

i	Magnetic field intensity Ii(nT)
		1	134.23
2	142.36
		3	157.68
4	148.19
		5	142.62
6	146.95
		7	162.27

Step 2: according to 7 positions of inertial navigation system instruction, reads respectively from the geomagnetic database prestored The geomagnetic field intensity information I (a of this position_i,b_i), and the gradient information I of the geomagnetic field intensity of this position_x,iAnd I_y,iSuch as table 3 institute Show:

Geomagnetic field intensity in table 3 geomagnetic database and gradient information

Step 3: initialization longitude and latitude error, course error and geomagnetic diurnal change error:

M=[0 00 0]^T

Step 4: calculate iterative parameter G, F and H according to formula (1), (2) and (3)

$G=\left[\begin{array}{c}-0.7557\\ -4.6567\\ 0.0495\\ -0.0116\end{array}\right]\×{10}^4$

$F=\left[\begin{array}{cccc}0.1878& 0.0494& 0.0237& 0.0008\\ 0.0494& 1.0446& -0.0328& 0.0014\\ 0.0237& -0.0328& 0.0092& 0.0001\\ 0.0008& 0.0014& 0.0001& 0.0001\end{array}\right]\×{10}^6$

$H=\left[\begin{array}{cccc}0.0161& 0.0013& 0.0201& 1.9347\\ 0.0013& 0.0023& 0.0112& 0.7456\\ 0.0201& 0.0112& 0.2207& 2.4284\\ 1.9347& 0.7456& 2.4284& 552.5414\end{array}\right]\×{10}^{-3}$

Step 5 is to step 9: choose ε=10^-6, default iterations is 50 times.Perform formula (4), formula (5), formula (7) to formula (10), iterative algorithm is implemented.And judge stopping criterion for iteration according to formula (6), it is known that when iterations is 7 Iteration ends.Longitude and latitude error, course error and geomagnetic diurnal change error M that iterative computation obtains be:

M=[-0.0325 0.0225-0.0382 16.2869]^T

Step 10: counted error amount substitutes into coupling track and reference locus relation equation (11), after calculating to mate Longitude, latitude value are as shown in table 4:

Table 4 matching result position coordinates

i	Longitude ui(°)	Latitude vi(°)
			1	127.1233	44.4205
2	127.1656	44.4395
			3	127.2275	44.4721
4	127.2882	44.5035
			5	127.3469	44.5344
6	127.4032	44.5657
			7	127.4556	44.5983

Above-mentioned iterative computation is the most time-consumingly 15.5 milliseconds.

Second iteration is used to improve the precision of matching locating method.During second iteration, by step one from inertial navigation system The position measurements a of the point to be matched read_iAnd b_iBy the coupling positioning result u obtained by first time iteration_iAnd v_iSubstitute, And geomagnetic diurnal change error initial in step 3 is substituted by the geomagnetic diurnal change error delta obtained by first time iteration.Second iteration Step is as follows:

Step one: substitute the point to be matched in 7 moment read from inertial navigation system with primary iteration result (being shown in Table 4) Position measurements a_iAnd b_i, as shown in table 5；Measured value according to proton magnetometer and magnetic field of the earth model, obtain current time And magnetic field intensity metrical information I in front 6 moment_i, as shown in table 2:

To be matched some position coordinates initial value during table 5 second iteration

i	Longitude ai(°)	Latitude bi(°)
			1	127.1233	44.4205
2	127.1656	44.4395
			3	127.2275	44.4721
4	127.2882	44.5035
			5	127.3469	44.5344
6	127.4032	44.5657
			7	127.4556	44.5983

Step 2: according to 7 positions in table 5, read these 7 positions from the geomagnetic database prestored respectively Geomagnetic field intensity I (a_i,b_i), and gradient I of the geomagnetic field intensity of this position_x,iAnd I_y,i, as shown in table 6:

Geomagnetic field intensity in table 6 geomagnetic database and gradient information

Step 3: initialize this geomagnetic diurnal change error by the geomagnetic diurnal change error delta obtained by for the first time iteration, and just Beginningization longitude and latitude error and course error:

M=[0 00 16.2869]^T

Step 4: calculate iterative parameter G, F and H according to formula (1), (2) and (3)

$G=\left[\begin{array}{c}-4.0307\\ 0.4355\\ -1.1574\\ -0.0353\end{array}\right]\×{10}^3$

$F=\left[\begin{array}{cccc}2.8352& -0.0371& 0.5985& 0.0102\\ -0.0371& 5.3856& -0.6101& -0.0016\\ 0.5985& -0.6101& 0.2257& 0.003\\ 0.0102& -0.0016& 0.003& 0.0001\end{array}\right]\×{10}^5$

$H=\left[\begin{array}{cccc}0.0174& -0.0078& -0.0709& 0.2682\\ -0.0078& 0.0076& 0.0540& -0.9626\\ -0.0709& 0.0540& 0.5067& -9.7977\\ 0.2682& -0.9626& -9.7977& 495.3714\end{array}\right]\×{10}^{-3}$

Step 5 is to step 9: choose ε=10^-6, default iterations is 50 times.Perform formula (4), formula (5), formula (7) to formula (10), iterative algorithm is implemented.And judge stopping criterion for iteration according to formula (6), it is known that when iterations is 7 Iteration ends.Longitude and latitude error, course error and geomagnetic diurnal change error M that iterative computation obtains be:

M=[0.0055-0.0086-0.0903 24.2326]^T

Step 10: counted error amount substitutes into coupling track and reference locus relation equation (11), after calculating to mate Longitude, latitude value are as shown in table 7:

Table 7 second iteration matching result position coordinates

i	Longitude ui(°)	Latitude vi(°)
			1	127.1288	44.4119
2	127.1693	44.4347
			3	127.2280	44.4727
4	127.2855	44.5094
			5	127.3413	44.5455
6	127.3945	44.5818
			7	127.4437	44.6190

Twice iterative computation is always time-consumingly 17.1 milliseconds in a computer.

In order to the effect after correcting geomagnetic diurnal change is described, we are not the most corrected matching result position during geomagnetic diurnal change Coordinate, as shown in table 8:

The matched position coordinate that table 8 corrects without geomagnetic diurnal change

i	Longitude ui′(°)	Latitude vi′(°)
			1	127.1841	44.4373
2	127.2299	44.4445
			3	127.2982	44.4596
4	127.3650	44.4739
			5	127.4298	44.4883
6	127.4923	44.5036
			7	127.5515	44.5213

For confirmatory experiment result, experiment sport car installs GPS navigation alignment system, such that it is able to obtain 7 moment True location coordinate, as shown in table 9:

Table 9 sport car true location coordinate

i	Longitude (°)	Latitude (°)
			1	127.1372	44.4094
2	127.1680	44.4314
			3	127.2232	44.4710
4	127.2786	44.5106
			5	127.3339	44.5502
6	127.3893	44.5898
			7	127.4447	44.6294

Positioning result (table 9) according to GPS navigation alignment system, can draw the location without diurnal correction in 7 moment Curve of error, there is the position error curve of the first time iteration result of diurnal correction and have the second time iteration result of diurnal correction Position error curve, as shown in Figures 2 and 3, wherein Fig. 2 gives longitude error curve, and it is bent that Fig. 3 gives latitude error Line.According to Fig. 2 and Fig. 3 it can be seen that the inertia earth magnetism matching locating method under the influence of the geomagnetic diurnal change proposed has higher Positioning precision, and after iteration, positioning precision can improve further.

The above, the only present invention preferably detailed description of the invention, these detailed description of the invention are all based on the present invention Different implementations under general idea, and protection scope of the present invention is not limited thereto, any are familiar with the art Technical staff in the technical scope that the invention discloses, the change that can readily occur in or replacement, all should contain the present invention's Within protection domain.Therefore, protection scope of the present invention should be as the criterion with the protection domain of claims.

Claims (1)

1. the inertia earth magnetism matching locating method under the influence of a geomagnetic diurnal change, it is characterised in that

Step one, read current time and the position measurements a of the point to be matched in front N-1 moment from inertial navigation system_iAnd b_i, Wherein a_iRepresent longitude, b_iRepresenting latitude, subscript i represents the most in the same time, i=1 ... N, N are integer and N ＞ 2, and i is that 1 expression is current In the moment, obtain current time and the measured value I of the geomagnetic field intensity in front N-1 moment by geomagnetic sensor_i；

Step 2, the position of the point N number of to be matched indicated according to inertial navigation system, read from the geomagnetic database prestored respectively Take the reference value I (a of corresponding geomagnetic field intensity_i,b_i), reference value I of the gradient of geomagnetic field intensity_x,iAnd I_y,i, wherein I_x,iTable Show the geomagnetic field intensity gradient in longitudinal value on i-th position, I_y,iRepresent that geomagnetic field intensity is latitudinal Gradient value on i-th position；

Step 3, introduce and initialize longitude and latitude error, course error and geomagnetic diurnal change error:

M=[Δ x Δ y α δ]^T=[0 00 0]^T

Step 4, calculate iterative parameter G, F and H according to formula (1), formula (2) and formula (3):

G=g (M) (1)

F=f (M) (2)

H=F^-1 (3)

Wherein:

G (M)=[g₁(M) g₂(M) g₃(M) g₄(M)]^T；

$f\left(M\right)=\left[\begin{array}{cccc}{f}_{11}\left(M\right)& f_{12}\left(M\right)& f_{13}\left(M\right)& f_{14}\left(M\right)\\ f_{12}\left(M\right)& f_{22}\left(M\right)& f_{23}\left(M\right)& f_{24}\left(M\right)\\ f_{13}\left(M\right)& f_{23}\left(M\right)& f_{33}\left(M\right)& f_{34}\left(M\right)\\ f_{14}\left(M\right)& f_{24}\left(M\right)& f_{34}\left(M\right)& f_{44}\left(M\right)\end{array}\right];$

$\begin{array}{c}{g}_1\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{x,i}\[I_{x,i}(a_i^{\′}\cos \mathrm{\α}+b_i^{\′}\sin \mathrm{\α}-a_i^{\′}+\mathrm{\Δ}x)\\ +I_{y,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α}-b_i^{\′}+\mathrm{\Δ}y)+I_{t,i}+\mathrm{\δ}\]\end{array};$

$\begin{array}{c}{g}_2\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{y,i}\[I_{x,i}(a_i^{\′}\cos \mathrm{\α}+b_i^{\′}\sin \mathrm{\α}-a_i^{\′}+\mathrm{\Δ}x)\\ +I_{y,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α}-b_i^{\′}+\mathrm{\Δ}y)+I_{t,i}+\mathrm{\δ}\]\end{array};$

$\begin{array}{c}{g}_3\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}\[I_{x,i}(a_i^{\′}\cos \mathrm{\α}+b_i^{\′}\sin \mathrm{\α}-a_i^{\′}+\mathrm{\Δ}x)+I_{y,i}(-a_i^{\′}\sin \mathrm{\α}\\ +b_i^{\′}\cos \mathrm{\α}-b_i^{\′}+\mathrm{\Δ}y)+I_{t,i}+\mathrm{\δ}\]\×\[I_{x,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α})\\ +I_{y,i}(-a_i^{\′}\cos \mathrm{\α}-b_i^{\′}\sin \mathrm{\α})\]\end{array};$

$\begin{array}{c}{g}_4\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}\[I_{x,i}\left(a_i^{\′}\cos \mathrm{\α}+b_i^{\′}\sin \mathrm{\α}-a_i^{\′}+\mathrm{\Δ}x\right)\\ +I_{y,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α}-b_i^{\′}+\mathrm{\Δ}y)+I_{t,i}+\mathrm{\δ}\]\end{array};$

$f_{11}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{x,i}^2;$

$f_{12}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{x,i}{I}_{y,i};$

$f_{13}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{x,i}\[I_{x,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α})+I_{y,i}(-a_i^{\′}\cos \mathrm{\α}-b_i^{\′}\sin \mathrm{\α})\];$

$f_{14}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{x,i};$

$f_{22}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{y,i}^2;$

$f_{23}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{y,i}\[I_{x,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α})+I_{y,i}(-a_i^{\′}\cos \mathrm{\α}-b_i^{\′}\sin \mathrm{\α})\];$

$f_{24}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{I}_{y,i};$

$\begin{array}{c}{f}_{33}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}{\[I_{x,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α})+I_{y,i}(-a_i^{\′}\cos \mathrm{\α}-b_i^{\′}\sin \mathrm{\α})\]}^2\\ +\[I_{x,i}(a_i^{\′}\cos \mathrm{\α}+b_i^{\′}\sin \mathrm{\α}-a_i^{\′}+\mathrm{\Δ}x)+I_{y,i}(-a_i^{\′}\sin \mathrm{\α}+b_i^{\′}\cos \mathrm{\α}\\ -b_i^{\′}+\mathrm{\Δ}y)+I_{t,i}+\mathrm{\δ}\]\×\[I_{x,i}(-a_i^{\′}\cos \mathrm{\α}-b_i^{\′}\sin \mathrm{\α})+I_{y,i}(a_i^{\′}\sin \mathrm{\α}\\ -b_i^{\′}\cos \mathrm{\α})\]\end{array};$

$f_{34}\left(M\right)=\underset{i=1}{\overset{N}{\Σ}}\[I_{x,i}(-a_i^{\′}sin\mathrm{\α}+b_i^{\′}cos\mathrm{\α})+I_{y,i}(-a_i^{\′}cos\mathrm{\α}-b_i^{\′}sin\mathrm{\α})\];$

f₄₄(M)=N；

a′_i=a_i-a₁, b '_i=b_i-b₁, I_t,i=I (a_i,b_i)-I_i；

Step 5, the calculating increment of longitude and latitude error, the increment of course error and the increment δ M of geomagnetic diurnal change error:

δ M=-H × G (4)

Step 6, renewal longitude and latitude error, course error and geomagnetic diurnal change error M:

M=M+ δ M (5)

Step 7, judge whether to meet and terminate iterated conditional, if meeting, stopping iteration and also jumping to step 10, otherwise jumping to step Eight；

Terminate iterated conditional for 1., 2. in any one or two: 1. iterations reaches preset times；2. longitude and latitude error 2 norms of the increment δ M of increment, the increment of course error and geomagnetic diurnal change error are less than setting value, i.e.

||δM||₂＜ ε (6)

Wherein ε is iteration minimum error set in advance；

Step 8, calculate parameter K and δ K according to the M after updating:

K=g (M) (7)

δ K=K-G (8)

Step 9, update iteration variable H and G according to formula (9) and formula (10), then jump to step 5,

H=H+ [(δ M)-H (δ K)] (δ M)^TH/[(δM)^TH(δK)] (9)

G=K (10)

Step 10, the longitude and latitude error that obtains according to iterative computation, course error and geomagnetic diurnal change error delta x, Δ y, α and δ, by institute Result substitutes into coupling track and i.e. obtains with the relation equation (11) of reference locus and mate track；

$\left[\begin{array}{c}{u}_i\\ v_i\end{array}\right]=\left[\begin{array}{cc}cos\mathrm{\α}& sin\mathrm{\α}\\ -sin\mathrm{\α}& cos\mathrm{\α}\end{array}\right]\left[\begin{array}{c}{a}_i-a_1\\ b_i-b_1\end{array}\right]+\left[\begin{array}{c}{a}_1\\ b_1\end{array}\right]+\left[\begin{array}{c}\mathrm{\Δ}x\\ \mathrm{\Δ}y\end{array}\right]---\left(11\right)$

Wherein u_iIt is the position longitude of the i-th moment matching result, v_iIt it is the position latitude of the i-th moment matching result.