CN105425304A - Compensation method for airplane aeromagnetic interference - Google Patents

Abstract

The invention provides a compensation method for airplane aeromagnetic interference, relates to the field of aeromagnetic detection, and aims at solving the problem that the estimation precision of the compensation coefficient is influenced by the morbid state of a present aeromagnetic interference coefficient calculation method. The method comprises the following steps that 1) a matrix delta is constructed; 2) the matrix delta is decomposed by wavelet; 3) bandbass filtering is carried out on a result of the step 2); 4) a subscript set L is determined so that the morbid state of bpf<L>(delta) is the weakest; and 5) a compensation coefficient equation is determined to solve the compensation coefficient of aeromagnetic interference. The method provided by the invention can effectively reduce the morbid state of the coefficient solution matrix, improve the solution precision of the compensation coefficient estimation algorithm, and can be applied to the field of aeromagnetic detection.

Description

A kind of compensation method of Aircraft magnetic interference

Technical field

The present invention relates to and belong to boat magnetic detection field, particularly relate to a kind of compensation method of Aircraft magnetic interference.

Background technology

Boat magnetic interference compensation technique is in boat magnetic detection process, remove a kind of technology of self magnetic interference of air transport platform.By analyzing type and the character of the magnetic interference of air transport platform self, set up air transport platform magnetic interference mathematical model, then in calibration flight process, record magnetic resultant field and three component seismic data according to the method for regulation, and be used for calculating the coefficient of air transport platform magnetic interference mathematical model.Actual carry out boat magnetic detection time, it also removes by the magnetic interference utilizing the model coefficient that estimates and aspect data estimation air transport platform to produce from magnetic resultant field, thus obtains not containing the magnetic field data of air transport platform magnetic interference.Mostly existing boat magnetic interference penalty coefficient computing method are based on T-L model, the magnetic interference of air transport platform is divided into stationary field, induction field and vortex field three types by this model, wherein stationary field coefficient has 3, induction field coefficient and vortex field coefficient respectively have 9, and both are relevant with the size of terrestrial magnetic field and rate of change respectively.Owing to the estimation procedure of penalty coefficient can be regarded as linear regression, therefore the character of coefficient solution matrix will certainly affect the precision of estimated result.Wherein, the pathosis of coefficient solution matrix is the key factor affecting penalty coefficient estimated accuracy.

Summary of the invention

In order to reduce the pathosis of coefficient solution matrix, the invention provides a kind of preprocess method reducing coefficient solution matrix pathosis, the method can be applicable to penalty coefficient algorithm for estimating and contributes to improving its solving precision.

This method solves the know-why that its technical matters adopts:

The magnetic field data that resultant field magnetometer records in calibration flight process is N × 1 column vector H _t, magnetic field value is N × 1 column vector H practically _e, and the magnetic interference that aircraft produces is N × 1 column vector H _i, then have:

H _T＝H _E+H _I(1)

Wherein, H is only had _tobtain by directly measuring;

The final purpose that boat magnetic interference compensates determines H _e, observe above formula and can find, as long as calculate H _iand by it from H _tin cut and just can obtain H _e.

According to T-L model, the magnetic interference that aircraft produces can be expressed as:

$H_I=\underset{i=1}{\overset{3}{\&Sigma;}}{p}_i{u}_i+\underset{i=1}{\overset{3}{\&Sigma;}}\underset{j=1}{\overset{3}{\&Sigma;}}{a}_{ij}{u}_i{u}_j+\underset{i=1}{\overset{3}{\&Sigma;}}\underset{j=1}{\overset{3}{\&Sigma;}}{b}_{ij}{\stackrel{\&CenterDot;}{u}}_i{u}_j---\left(2\right)$

Wherein, p _i, a _ij, b _ijfor needing the boat magnetic interference penalty coefficient estimated, u _iand u _jfor the direction cosine by terrestrial magnetic field and the axial angle of aircraft, for u _iderivative;

The magnetic field three-component that three-component magnetometer records is v ₁, v ₂, v ₃, u _ican calculate according to following formula:

$u_i=\frac{v_i}{\sqrt{v_1^2+v_2^2+v_3^2}},(i=1,2,3)---\left(3\right)$

According to geometrical symmetry, formula (2) can change into more compact form:

$\begin{array}{c}{H}_I=p_1{u}_1+p_2{u}_2+p_3{u}_3+a_{11}{u}_1^2+a_{12}{u}_1{u}_2+a_{13}{u}_1{u}_3+a_{22}{u}_2^2+a_{23}{u}_2{u}_3\\ +b_{11}{\stackrel{\&CenterDot;}{u}}_1{u}_1+b_{12}{\stackrel{\&CenterDot;}{u}}_1{u}_2+b_{13}{\stackrel{\&CenterDot;}{u}}_1{u}_3+b_{21}{\stackrel{\&CenterDot;}{u}}_2{u}_1+b_{22}{\stackrel{\&CenterDot;}{u}}_2{u}_2+b_{23}{\stackrel{\&CenterDot;}{u}}_2{u}_3+b_{31}{\stackrel{\&CenterDot;}{u}}_3{u}_1+b_{32}{\stackrel{\&CenterDot;}{u}}_3{u}_3\\ =\mathrm{\&delta;}\mathrm{\&theta;}\end{array}---\left(4\right)$

In formula, θ is the column vector be made up of undetermined coefficient, and δ is by u _i, u _iu _jwith the row vector formed.

According to formula (4), H _ican be expressed as:

H _I＝Δθ(5)

Wherein, for N × 16 matrix; The magnetic interference that aircraft produces is relevant with the frequency of calibration maneuver, and the frequency of terrestrial magnetic field and calibration maneuver has nothing to do, if bpf is FIR bandpass filter, convolution (1) and (2) then have:

bpf(H _T)＝bpf(H _E)+bpf(H _I)＝bpf(H _E)+bpf(Δ)θ(6)

Wherein bpf () expression carries out bandpass filtering to each row of matrix;

Free transmission range due to bpf calibrates the set of frequency of maneuver according to aircraft, therefore bpf (H _e) ≈ 0, then formula (6) is converted into:

bpf(H _I)＝bpf(Δ)θ(7)

The column vector that penalty coefficient is formed can calculate according to formula (7).

Easily find out the estimation that the character of bpf (Δ) will directly affect θ.Because Δ knows, the character of bpf (Δ) therefore can be affected by the correlation parameter of adjustment bpf.

It is a principal object of the present invention to the pathosis degree reducing bpf (Δ), utilize Wavelet Decomposition Technology to make the pathosis of bpf (Δ) minimum;

The present invention solves because the pathosis of coefficient solution matrix affects the problem of penalty coefficient estimated accuracy in existing boat magnetic interference penalty coefficient computing method, and proposes a kind of compensation method of Aircraft magnetic interference.

A compensation method for Aircraft magnetic interference, carry out according to the following steps:

One, structural matrix Δ:

The magnetic field data that resultant field magnetometer records in calibration flight process is N × 1 column vector H _t, magnetic field value is N × 1 column vector H practically _e, and the magnetic interference that aircraft produces is N × 1 column vector H _i, then have:

H _T＝H _E+H _I(1)

Wherein, H is only had _tobtain by directly measuring;

The final purpose that boat magnetic interference compensates determines H _e, observe above formula and can find, as long as calculate H _iand by it from H _tin cut and just can obtain H _e;

According to T-L model, the magnetic interference that aircraft produces can be expressed as:

$H_I=\underset{i=1}{\overset{3}{\&Sigma;}}{p}_i{u}_i+\underset{i=1}{\overset{3}{\&Sigma;}}\underset{j=1}{\overset{3}{\&Sigma;}}{a}_{ij}{u}_i{u}_j+\underset{i=1}{\overset{3}{\&Sigma;}}\underset{j=1}{\overset{3}{\&Sigma;}}{b}_{ij}{\stackrel{\&CenterDot;}{u}}_i{u}_j---\left(2\right)$

Wherein, p _i, a _ij, b _ijfor needing the boat magnetic interference penalty coefficient estimated, u _iand u _jfor the direction cosine by terrestrial magnetic field and the axial angle of aircraft, for u _iderivative;

The magnetic field three-component that three-component magnetometer records is v ₁, v ₂, v ₃, u _ican calculate according to following formula:

$u_i=\frac{v_i}{\sqrt{v_1^2+v_2^2+v_3^2}},(i=1,2,3)---\left(3\right)$

According to geometrical symmetry, formula (2) can change into more compact form:

$\begin{array}{c}{H}_I=p_1{u}_1+p_2{u}_2+p_3{u}_3+a_{11}{u}_1^2+a_{12}{u}_1{u}_2+a_{13}{u}_1{u}_3+a_{22}{u}_2^2+a_{23}{u}_2{u}_3\\ +b_{11}{\stackrel{\&CenterDot;}{u}}_1{u}_1+b_{12}{\stackrel{\&CenterDot;}{u}}_1{u}_2+b_{13}{\stackrel{\&CenterDot;}{u}}_1{u}_3+b_{21}{\stackrel{\&CenterDot;}{u}}_2{u}_1+b_{22}{\stackrel{\&CenterDot;}{u}}_2{u}_2+b_{23}{\stackrel{\&CenterDot;}{u}}_2{u}_3+b_{31}{\stackrel{\&CenterDot;}{u}}_3{u}_1+b_{32}{\stackrel{\&CenterDot;}{u}}_3{u}_3\\ =\mathrm{\&delta;}\mathrm{\&theta;}\end{array}---\left(4\right)$

In formula, θ is the column vector be made up of undetermined coefficient, and δ is by u _i, u _iu _jwith the row vector formed;

According to formula (4), H _ican be expressed as:

H _I＝Δθ(5)

Wherein, $\mathrm{\&Delta;}={({\mathrm{\&delta;}}_1^T,{\mathrm{\&delta;}}_2^T,...,{\mathrm{\&delta;}}_N^T)}^T$ For N × 16 matrix;

Two, wavelet decomposition matrix Δ is utilized:

Each of matrix Δ is classified as f _i, (i=1,2 ... 16), then utilize wavelet decomposition can by f _iresolve into following form:

f _i＝f ₀+w ₁+w ₂+…+w _j(8)

W in above formula _l, (l=1,2 ..., j) represent f _iin different frequency composition, and the larger w of l _lfrequency higher, w _jfrequency the highest expression f _iin noise; f ₀represent f _iin extremely low frequency composition, i.e. geomagnetic field component;

Three, bandpass filtering is realized to the result of step 2:

Definition U _i={ w _il| 1≤l≤j-1} is the result of the i-th row in matrix Δ being carried out to wavelet decomposition, wherein f _i0and w _ijall ignore; For arbitrary U _ineed to determine identical continuous indexed set L={s, s+1 ... t} (1≤s≤t≤j-1), then for f _irespective belt pass filtering results be:

bpf _L(f _i)＝w _is+w _i(s+1)+…w _it(9)

And then have:

bpf _L(Δ)＝(bpf _L(f ₁),bpf _L(f ₂),…,bpf _L(f ₁₆))(10)

Four, by determining that indexed set L makes bpf _lthe pathosis of (Δ) is the most weak:

Conditional number is utilized to carry out metric matrix bpf _lthe pathosis degree of (Δ):

$cond\left({\mathrm{bpf}}_L\right(\mathrm{\&Delta;}\left)\right)=\frac{{\mathrm{\&lambda;}}_{max}}{{\mathrm{\&lambda;}}_{min}}---\left(11\right)$

Wherein, λ _maxand λ _minrepresent bpf respectively _lthe eigenvalue of maximum of (Δ) and minimal eigenvalue, cond (bpf _l(Δ)) value is larger, and the pathosis of matrix is more serious;

Owing to needing the subscript chosen to be continuous print, known L={s, s+1 ... total (j-1) j/2 kind different situations of value of t} (1≤s≤t≤j-1), the value of j is no more than 10, utilizes the method for exhaustion to determine L in the present invention _opt, make minimum.

Five, determine penalty coefficient equation and solve boat magnetic interference penalty coefficient:

Coefficient equation is compensated as follows according to step 4:

${\mathrm{bpf}}_{L_{opt}}\left(H_I\right)={\mathrm{bpf}}_{L_{opt}}\left(\mathrm{\&Delta;}\right)\mathrm{\&theta;}---\left(12\right)$

Boat magnetic interference penalty coefficient is solved by above formula.

The present invention includes following beneficial effect:

1, the method proposed effectively reduces the pathosis of coefficient solution matrix, and then improves the solving precision of penalty coefficient algorithm for estimating;

2, the inventive method, in order to realize bandpass filtering, reasonably have ignored w _jcomponent and f ₀component, neither affects result of calculation, in turn simplify computation process;

3, utilize the method for exhaustion to determine L in the present invention _opt, make minimum, farthest reduce the pathosis of coefficient solution matrix.

Embodiment

For enabling above-mentioned purpose of the present invention, feature and advantage become apparent more, and below in conjunction with embodiment, the present invention is further detailed explanation.

The compensation method of a kind of Aircraft magnetic interference described in embodiment one, present embodiment, carry out according to the following steps:

One, structural matrix Δ:

The magnetic field data that resultant field magnetometer records in calibration flight process is N × 1 column vector H _t, magnetic field value is N × 1 column vector H practically _e, and the magnetic interference that aircraft produces is N × 1 column vector H _i, then have:

H _T＝H _E+H _I(1)

Wherein, H is only had _tobtain by directly measuring;

The final purpose that boat magnetic interference compensates determines H _e, observe above formula and can find, as long as calculate H _iand by it from H _tin cut and just can obtain H _e;

According to T-L model, the magnetic interference that aircraft produces can be expressed as:

$H_I=\underset{i=1}{\overset{3}{\&Sigma;}}{p}_i{u}_i+\underset{i=1}{\overset{3}{\&Sigma;}}\underset{j=1}{\overset{3}{\&Sigma;}}{a}_{ij}{u}_i{u}_j+\underset{i=1}{\overset{3}{\&Sigma;}}\underset{j=1}{\overset{3}{\&Sigma;}}{b}_{ij}{\stackrel{\&CenterDot;}{u}}_i{u}_j---\left(2\right)$

Wherein, p _i, a _ij, b _ijfor needing the boat magnetic interference penalty coefficient estimated, u _iand u _jfor the direction cosine by terrestrial magnetic field and the axial angle of aircraft, for u _iderivative;

The magnetic field three-component that three-component magnetometer records is v ₁, v ₂, v ₃, u _ican calculate according to following formula:

$u_i=\frac{v_i}{\sqrt{v_1^2+v_2^2+v_3^2}},(i=1,2,3)---\left(3\right)$

According to geometrical symmetry, formula (2) can change into more compact form:

$\begin{array}{c}{H}_I=p_1{u}_1+p_2{u}_2+p_3{u}_3+a_{11}{u}_1^2+a_{12}{u}_1{u}_2+a_{13}{u}_1{u}_3+a_{22}{u}_2^2+a_{23}{u}_2{u}_3\\ +b_{11}{\stackrel{\&CenterDot;}{u}}_1{u}_1+b_{12}{\stackrel{\&CenterDot;}{u}}_1{u}_2+b_{13}{\stackrel{\&CenterDot;}{u}}_1{u}_3+b_{21}{\stackrel{\&CenterDot;}{u}}_2{u}_1+b_{22}{\stackrel{\&CenterDot;}{u}}_2{u}_2+b_{23}{\stackrel{\&CenterDot;}{u}}_2{u}_3+b_{31}{\stackrel{\&CenterDot;}{u}}_3{u}_1+b_{32}{\stackrel{\&CenterDot;}{u}}_3{u}_3\\ =\mathrm{\&delta;}\mathrm{\&theta;}\end{array}---\left(4\right)$

In formula, θ is the column vector be made up of undetermined coefficient, and δ is by u _i, u _iu _jwith the row vector formed;

According to formula (4), H _ican be expressed as:

H _I＝Δθ(5)

Wherein, for N × 16 matrix;

Two, wavelet decomposition matrix Δ is utilized:

Utilize Daubechies small echo to each row f of Δ _icarry out wavelet decomposition, the number of plies of decomposition is 10, and namely the value of j is 10, f _iresolve into following form:

f _i＝f ₀+w ₁+w ₂+…+w ₁₀

W in above formula _l, (l=1,2 ..., 10) and represent f _iin different frequency composition, and the larger w of l _lfrequency higher, w ₁₀frequency the highest expression f _iin noise; f ₀represent f _iin extremely low frequency composition, i.e. geomagnetic field component;

Three, bandpass filtering is realized to the result of step 2:

Definition U _i={ w _il| 1≤l≤9} is the result of the i-th row in matrix Δ being carried out to wavelet decomposition, wherein f _i0and w _i10all ignore;

For arbitrary U _ineed to determine identical continuous indexed set L={s, s+1 ... t} (1≤s≤t≤9), then for f _irespective belt pass filtering results be:

bpf _L(f _i)＝w _is+w _i(s+1)+…w _it

And then have:

bpf _L(Δ)＝(bpf _L(f ₁),bpf _L(f ₂),…,bpf _L(f ₁₆))

Four, by determining that indexed set L makes bpf _lthe pathosis of (Δ) is the most weak:

Conditional number is utilized to carry out metric matrix bpf _lthe pathosis degree of (Δ):

$cond\left({\mathrm{bpf}}_L\right(\mathrm{\&Delta;}\left)\right)=\frac{{\mathrm{\&lambda;}}_{max}}{{\mathrm{\&lambda;}}_{min}}$

Wherein, λ _maxand λ _minrepresent bpf respectively _lthe eigenvalue of maximum of (Δ) and minimal eigenvalue, cond (bpf _l(Δ)) value is larger, and the pathosis of matrix is more serious;

Owing to needing the subscript chosen to be continuous print, known L={s, s+1 ... the value of t} (1≤s≤t≤9) has 9 × 10/2 kinds of different situations, utilizes the method for exhaustion to determine L in the present invention _opt, make minimum.

Five, determine penalty coefficient equation and solve boat magnetic interference penalty coefficient:

Coefficient equation is compensated as follows according to step 4:

${\mathrm{bpf}}_{L_{opt}}\left(H_I\right)={\mathrm{bpf}}_{L_{opt}}\left(\mathrm{\&Delta;}\right)\mathrm{\&theta;}$

Boat magnetic interference penalty coefficient is solved by above formula.

Embodiment two, present embodiment and embodiment one unlike: the number of plies of decomposing described in step 2 is 9, and namely the value of j is 9.Other step is identical with embodiment one with parameter.

Embodiment three, present embodiment and embodiment one unlike: the number of plies of decomposing described in step 2 is 8, and namely the value of j is 8.Other step is identical with embodiment one with parameter.

Embodiment four, present embodiment and embodiment one unlike: the number of plies of decomposing described in step 2 is 7, and namely the value of j is 7.Other step is identical with embodiment one with parameter.

Claims (4)

1. a compensation method for Aircraft magnetic interference, is characterized in that it carries out according to the following steps:

One, structural matrix Δ:

The magnetic field data that resultant field magnetometer records in calibration flight process is N × 1 column vector H _t, magnetic field value is N × 1 column vector H practically _e, and the magnetic interference that aircraft produces is N × 1 column vector H _i, then have:

H _T＝H _E+H _L(1)

Wherein, H _tobtain by directly measuring;

The final purpose that boat magnetic interference compensates determines H _e, observe above formula and can find, as long as calculate H _iand by it from H _tin cut and just can obtain H _e;

According to T-L model, the magnetic interference that aircraft produces can be expressed as:

$H_I=\underset{i=1}{\overset{3}{\&Sigma;}}{p}_i{u}_i+\underset{i=1}{\overset{3}{\&Sigma;}}\underset{j=1}{\overset{3}{\&Sigma;}}{a}_{ij}{u}_i{u}_j+\underset{i=1}{\overset{3}{\&Sigma;}}\underset{j=1}{\overset{3}{\&Sigma;}}{b}_{ij}{\stackrel{\&CenterDot;}{u}}_i{u}_j---\left(2\right)$

Wherein, p _i, a _ij, b _ijfor needing the boat magnetic interference penalty coefficient estimated, u _iand u _jfor the direction cosine by terrestrial magnetic field and the axial angle of aircraft, for u _iderivative;

The magnetic field three-component that three-component magnetometer records is v ₁, v ₂, v ₃, u _ican calculate according to following formula:

$u_i=\frac{v_i}{\sqrt{v_1^2+v_2^2+v_3^2}},(i=1,2,3)---\left(3\right)$

According to geometrical symmetry, formula (2) can be converted into following formula:

$\begin{array}{c}{H}_I=p_1{u}_1+p_2{u}_2+p_3{u}_3+a_{11}{u}_1^2+a_{12}{u}_1{u}_2+a_{13}{u}_1{u}_3+a_{22}{u}_2^2+a_{23}{u}_2{u}_3\\ +b_{11}{\stackrel{\&CenterDot;}{u}}_1{u}_1+b_{12}{\stackrel{\&CenterDot;}{u}}_1{u}_2+b_{13}{\stackrel{\&CenterDot;}{u}}_1{u}_3+b_{21}{\stackrel{\&CenterDot;}{u}}_2{u}_1+b_{22}{\stackrel{\&CenterDot;}{u}}_2{u}_2+b_{23}{\stackrel{\&CenterDot;}{u}}_2{u}_3+b_{31}{\stackrel{\&CenterDot;}{u}}_3{u}_1+b_{32}{\stackrel{\&CenterDot;}{u}}_3{u}_3\\ =\mathrm{\&delta;}\mathrm{\&theta;}\end{array}---\left(4\right)$

In formula, θ is the column vector be made up of undetermined coefficient, and δ is by u _i, u _iu _jwith the row vector formed;

According to formula (4), H _ican be expressed as:

H _I＝Δθ(5)

Wherein, $\mathrm{\&Delta;}={({\mathrm{\&delta;}}_1^T,{\mathrm{\&delta;}}_2^T,...,{\mathrm{\&delta;}}_N^T)}^T$ For N × 16 matrix;

Two, wavelet decomposition matrix Δ is utilized:

Each of matrix Δ is classified as f _i, (i=1,2 ... 16), then utilize wavelet decomposition can by f _iresolve into following form:

f _i＝f ₀+w ₁+w ₂+…+w _j(8)

W in above formula _l, (l=1,2 ..., j) represent f _iin different frequency composition, and the larger w of l _lfrequency higher, w _jfrequency the highest expression f _iin noise; f ₀represent f _iin extremely low frequency composition, i.e. geomagnetic field component;

Three, bandpass filtering is realized to the result of step 2:

Definition U _i={ w _il| 1≤l≤j-1} is the result of the i-th row in matrix Δ being carried out to wavelet decomposition, wherein f _i0and w _ijall ignore; For arbitrary U _ineed to determine identical continuous indexed set L={s, s+1 ... t} (1≤s≤t≤j-1), then for f _irespective belt pass filtering results be:

bpf _L(f _i)＝w _is+w _i(s+1)+…w _it(9)

And then have:

bpf _L(Δ)＝(bpf _L(f ₁),bpf _L(f ₂),…,bpf _L(f ₁₆))(10)

Four, by determining that indexed set L makes bpf _lthe pathosis of (Δ) is the most weak:

Conditional number is utilized to carry out metric matrix bpf _lthe pathosis degree of (Δ):

$\mathrm{cond}\left({\mathrm{bpf}}_L\left(\mathrm{\&Delta;}\right)\right)=\frac{{\mathrm{\&lambda;}}_{\max }}{{\mathrm{\&lambda;}}_{\min }}---\left(11\right)$

Wherein, λ _maxand λ _minrepresent bpf respectively _lthe eigenvalue of maximum of (Δ) and minimal eigenvalue, cond (bpf _l(Δ)) value is larger, and the pathosis of matrix is more serious;

Five, determine penalty coefficient equation and solve boat magnetic interference penalty coefficient:

Coefficient equation is compensated as follows according to step 4:

${\mathrm{bpf}}_{L_{\mathrm{opt}}}\left(H_I\right)={\mathrm{bpf}}_{L_{\mathrm{opt}}}\left(\mathrm{\&Delta;}\right)\mathrm{\&theta;}---\left(12\right)$

Wherein, L in indexed set _optmake value is minimum, solves boat magnetic interference penalty coefficient by above formula.

2. the compensation method of a kind of Aircraft magnetic interference as claimed in claim 1, it is characterized in that the number of plies of decomposing described in step 2 is 9, namely the value of j is 9.

3. the compensation method of a kind of Aircraft magnetic interference as claimed in claim 1, it is characterized in that the number of plies of decomposing described in step 2 is 8, namely the value of j is 8.

4. the compensation method of a kind of Aircraft magnetic interference as claimed in claim 1, it is characterized in that the number of plies of decomposing described in step 2 is 7, namely the value of j is 7.