CN103792550B - A kind of associating anti-interference method based on array antenna and GPS/SINS - Google Patents

Abstract

The present invention is to provide a kind of associating anti-interference method based on array antenna and GPS/SINS.After the position of initialization carrier and attitude, set up GPS/SINS integrated navigation state equation and measurement equation; GPS/SINS integrated navigation provides position and the attitude of carrier in real time, calculates the position of present satellites according to satellite ephemeris information, obtains the steering vector between satellite to carrier; Described steering vector as the prior imformation of multiple constraint minimum variance space-time adaptive Processing Algorithm, in spatial domain, time domain suppress simultaneously broadband interference and arrowband interference.The present invention can form wave beam in many satellites in view directions simultaneously, forms zero and falls into, thus suppress undesired signal while strengthening satellite-signal at interference radiating way.The present invention adopts circular configuration aerial array, and the Wave beam forming that GPS/SINS integrated navigation is array antenna provides position and the attitude of carrier, adopts satellite ephemeris to provide the position of satellite, thus provides prior imformation for Wave beam forming.

Description

A kind of associating anti-interference method based on array antenna and GPS/SINS

Technical field

What the present invention relates to is a kind of integrated navigation technology, specifically a kind of array antenna and GPS(GPS)/SINS(strapdown inertial navigation system) integrated navigation technology.

Background technology

GPS navigation system round-the-clockly can provide exact position covering the whole world, speed and temporal information for its user, and its using value is more and more higher.But along with the raising of artificial interference technology, satellite only relies on its spread spectrum system to carry out anti-interferencely can not meeting consumers' demand.According to ICD-200, the anti-jamming margin of commercial GPS receivers is no more than 24dB(and depends on noise level).If namely jamming-to-signal ratio is greater than 24dB, commercial GPSC/A code receiver just cannot keep the tracking to signal.Test shows, power is that the jammer of 1W can make the C/A code receiver within 85 kilometers work.In addition, existing gps signal is at well-known transmitted on frequencies, and its modulation signature is widely known by the people, and signal to noise ratio (S/N ratio) is lower again, is thus easy to carry out disturbing or cheating.Therefore, study GPS Anti-Jamming Technique and become focus and emphasis.

At present, the anti-interference method based on signal transacting studies field the most active, and common technology comprises temporal filtering technique and airspace filter technology.Adaptive temporal filter technology is a kind of narrow-band interference rejection method, and the useful signal received, interference and noise make required cost function minimize by adaptive algorithm and remove undesired signal by it.Time-domain filtering is for being very stable during BOUNDED DISTURBANCES source, this is because it can provide complicated rejection filter criterion simultaneously, and it is counted as the independent embedded part before the embedded part of GPS pre-process and post-process or GPS.When very little to size impact, this technology is greater than 30dB to the suppression that arrowband disturbs.For the interference of elimination arrowband, temporal filtering technique can be used for complicated arrowband and continuous wave CO_2 laser source, but this can be subject to again the restriction of remaining computation bandwidth, and this remaining computation bandwidth can hinder effective gps signal process.Airspace filter technology is realized by adaptive array, and it can effectively suppress coherent interference and broadband interference.Adaptive antenna and self-adaptive filters in time area somewhat similar, also certain cost function will be made to minimize, if it have a very large defect be exactly the incident direction of useful signal and undesired signal close to each other time, adaptive antenna also can impact useful signal while suppression undesired signal.Its realization needs to consume larger power, more spend, and the workbench that needs one are larger.Airspace filter mainly contains two types: zero falls into and Wave beam forming.The demand of null steering technique to useful signal quantity of information is minimum, and it is easier than beam-forming technology realizes.The output power of null steering technique minimum signal, it can be divided into arrowband and two kinds, broadband situation.Beam-forming technology requires the information of many relevant useful signals, and implements also more complicated.The output signal-to-noise ratio of beam-forming technology maximum signal, beam-forming technology is also divided into arrowband and two kinds, broadband situation.

Simple time-domain filtering and simple spatial domain have relative merits separately, but the quality of these two kinds of disposal routes just can be complementary, by the two connected applications, and then a kind of associating Anti-Jamming Technique, i.e. space-time two-dimension Combined Treatment (STAP) Anti-Jamming Technique can be formed.The airspace filter that the spatial processing energy force rate of space-time two-dimensional Combined Treatment is simple is stronger.When needs to broadband signal carry out zero fall into time, because signal bandwidth be can not ignore, what is called can be caused in space " to disperse " phenomenon.General airspace filter means all need the space of causing this broadband " to disperse " to consider especially, and this is also the reason that simple spatial domain method process broadband signal is difficult to obtain promising result.Frequency information and spatial information (si) are combined by space-time two-dimensional Combined Treatment, complete the showing of each component of signal naturally when sky, plane will be able to dispersed, thus can eliminate interference to greatest extent and suppress.Further, STAP technology also has the plurality of advantages such as inherent Wave beam forming, the signal equalization of inherence and the anti-multipath jamming ability of inherence, thus can realize anti-interference process while enhancing signal.

In addition, GPS/SINS integrated navigation also can improve the antijamming capability of GPS.The loop bandwidth of receiver needs tradeoff design between antijamming capability and performance of dynamic tracking, namely receiver is in order to improve the ability that self suppresses interference and noise, need loop bandwidth to reduce, and need loop bandwidth to increase in order to the dynamic property of following the tracks of carrier.And after introducing Inertia information, SINS accurately can estimate the speed of carrier, calculate the Doppler shift of carrier relative to satellite, thus reduce the loop bandwidth of receiver, increase the antijamming capability of receiver.

But all there is certain defect in above traditional anti-interference method.Simple time-domain processing method can not suppress broadband interference, and zero in airspace filter falls into method can not provide gain in satellite-signal direction, and satellite-signal is also likely suppressed.Beamforming Method in airspace filter and space-time adaptive process often need satellite to arrive the angle of receiver.GPS/SINS combines anti-interference method needs GPS correctly can export pseudorange, pseudorange rates or position, velocity information.

Summary of the invention

The object of the present invention is to provide a kind of associating anti-interference method based on array antenna and GPS/SINS suppressing interference performance strong.

The object of the present invention is achieved like this:

After the position of initialization carrier and attitude, set up GPS/SINS integrated navigation state equation and measurement equation; GPS/SINS integrated navigation provides position and the attitude of carrier in real time, calculates the position of present satellites according to satellite ephemeris information simultaneously, thus obtains the steering vector between satellite to carrier, and namely satellite arrives deflection and the angle of pitch of carrier; Then, described steering vector as the prior imformation of multiple constraint minimum variance space-time adaptive process (MCMV-STAP) algorithm, in spatial domain, time domain suppress simultaneously broadband interference and arrowband interference.Specifically comprise the steps:

Step 1: initialization carrier positions, speed and attitude information

In earth coordinates, the coordinate of setting carrier initial time: latitude L, longitude λ and height h; The speed of initialization carrier northeastward in sky coordinate system: east orientation speed V _e, north orientation speed V _nwith sky to speed V _u; The attitude angle of initialization carrier, comprises pitching angle theta, roll angle γ and position angle ψ; Then set the flight path of carrier, static, the track such as rectilinear motion or circular motion can be set as.Thus obtain the output of desirable gyroscope and accelerometer, f _e, f _nand f _ube expressed as accelerometer in east orientation, north orientation and sky to output.

Step 2: the parameter information of initialization Navigation Filter

In GPS/SINS integrated navigation wave filter, adopt feedback compensation mode.The quantity of state of Navigation Filter is attitude error, velocity error, site error, gyro error and accelerometer error.φ _e, φ _nand φ _ube expressed as carrier angle of pitch error, roll angle error and azimuth angle error, δ V _e, δ V _nwith δ V _ube expressed as the east orientation velocity error of carrier, north orientation velocity error and sky to velocity error, δ L, δ λ and δ h are expressed as the latitude error of carrier, longitude error and height error.ε _e, ε _nand ε _ube expressed as gyroscope in east orientation, north orientation and sky to drift. with be expressed as accelerometer in east orientation, north orientation and sky to output error.

Step 3: the error rate of computing gyroscope and accelerometer

Step 4: according to the parameter in step 1-3, calculates attitude of carrier angle error rate of change, velocity error rate of change and site error rate of change

Step 5: introduce GPS pseudo range measurement information, adopts feedback compensation mode to correct SINS output information, obtains current position accurately and attitude.

Code phase error measured by GPS track loop calculates the pseudorange between satellite to carrier, and pseudorange upgrades the quantity of state of Navigation Filter as measurement information, thus dopes the error of current SINS institute state quantity measurement amount.Then adopt feedback compensation mode to correct SINS output information, obtain current position accurately and attitude.

Step 6: calculate the coordinate of carrier in ECEF coordinate system

The carrier positions calculated by step 5 is arranged in earth coordinates, i.e. (L, λ, h), is translated in ECEF coordinate system, its position coordinates can be expressed as

$X_p^e={[(R+h)\cos L\cos \mathrm{\λ},(R+h)\cos L\sin \mathrm{\λ},(R+h)\sin L]}^T---\left(1\right)$

In formula (1), R is earth radius.

Step 7: calculate carrier coordinate system (b system) is tied to ECEF coordinate system (e system) transition matrix to the transition matrix of navigational coordinate system (n system), navigation coordinate

The attitude of carrier angle calculated by step 5, can calculate the transition matrix of carrier coordinate system to navigational coordinate system for

$C_b^n=\left[\begin{array}{ccc}\cos \mathrm{\γ}\cos \mathrm{\ψ}-\sin \mathrm{\ψ}\sin \mathrm{\θ}\sin \mathrm{\γ}& -\sin \mathrm{\ψ}\cos \mathrm{\θ}& \sin \mathrm{\γ}\cos \mathrm{\ψ}+\cos \mathrm{\γ}\sin \mathrm{\θ}\sin \mathrm{\ψ}\\ \cos \mathrm{\γ}\sin \mathrm{\ψ}+\sin \mathrm{\γ}\sin \mathrm{\θ}\cos \mathrm{\ψ}& \cos \mathrm{\ψ}\cos \mathrm{\θ}& \sin \mathrm{\γ}\sin \mathrm{\ψ}-\cos \mathrm{\γ}\sin \mathrm{\θ}\cos \mathrm{\ψ}\\ -\sin \mathrm{\γ}\cos \mathrm{\θ}& \sin \mathrm{\θ}& \cos \mathrm{\θ\γ}\cos \end{array}\right]---\left(2\right)$

The carrier positions calculated by step 5, can calculate the transition matrix that navigation coordinate is tied to ECEF coordinate system for

$C_n^e=\left[\begin{array}{ccc}-\sin \mathrm{\λ}& -\sin L\cos \mathrm{\λ}& \cos L\cos \mathrm{\λ}\\ \cos \mathrm{\λ}& -\sin L\sin \mathrm{\λ}& \cos L\sin \mathrm{\λ}\\ 0& \cos L& \sin L\end{array}\right]---\left(3\right)$

Step 8: calculate the steering vector between carrier to satellite

Can in the hope of the transformed matrix of carrier coordinate system to ECEF coordinate system by (2), (3) in step 7

$C_b^e=C_n^e{C}_b^n---\left(4\right)$

Satellite ephemeris is utilized to calculate the position of satellite in ECEF coordinate system and in conjunction with formula (1) (4), can in the hope of the steering vector between carrier to satellite

${\stackrel{\→}{r}}^b={\left(C_b^e\right)}^T(X_s^e-X_p^e)---\left(5\right)$

Step 9: in carrier coordinate system, calculates position angle and the angle of pitch that satellite arrives antenna

If will coordinate definition in carrier coordinate system is then satellite arrives the azimuth angle alpha of antenna and angle of pitch β can be expressed as respectively

$\mathrm{\α}=\arctan (\hat{x},\hat{y})---\left(6\right)$

$\mathrm{\β}=\arctan (\hat{z},\sqrt{{\hat{x}}^2+{\hat{y}}^2})---\left(7\right)$

Step 10: design circular array antenna structure

In order to control beam position satellites in view direction at position angle and angle of pitch direction simultaneously, the present invention adopts circular array antenna structure.6 array element is uniformly distributed in round battle array circumferentially, and make radius of circle be r, then the interval circumferentially between adjacent array element is also r.Choosing of array element interval, and to meet Nyquist's theorem the same at time-domain sampling interval, and Space domain sampling interval d should be less than 1/2 of satellite carrier wavelength X.By satellite frequency f=1575.42 × 10 ⁶mHz, so array element distance is

$d\≤\frac{\mathrm{\λ}}{2}=\frac{1}{2}\·\frac{c}{f_{l1}}=\frac{1}{2}\·\frac{3\×{10}^8}{1575.42\×{10}^6}=\frac{1}{2}\·0.19=0.095m=9.5\mathrm{cm}---\left(8\right)$

In formula (8), c is the light velocity.

In order to make main beam width narrower, secondary lobe is lower, and resolution is high, and what need to be tried one's best in array element interval is large, so get radius of circle r=d=9.5cm, the diameter of whole array antenna is about 19cm.

Step 11: calculate the time delay between each array element of satellite-signal arrival antenna

According to (6) (7) in step 9, can by α and β representation unit vector

e(α,β)=(sinαcosβ,cosαcosβ,sinβ) ^T(9)

Therefore, satellite-signal arrives i-th bay and to the mistiming τ between first reference array element _ican be expressed as

τ _i=e ^T·(x _i-x ₁)/ci=1,2,…M-1(10)

In formula (10), M representative antennas array element number.

Step 12: set up array antenna received signals model

User generally can receive the satellite-signal of more than 4, therefore, needs to form satellite corresponding to multiple beam position.Suppose that array antenna received has arrived P satellite-signal, Q undesired signal, then the signal model that antenna receives can be expressed as

$U\left(t\right)=\underset{k=1}{\overset{P}{\mathrm{\Σ}}}a({\mathrm{\α}}_k,{\mathrm{\β}}_k,T_s)s_k\left(t\right)+\underset{k=P+1}{\overset{P+Q}{\mathrm{\Σ}}}a({\mathrm{\α}}_k,{\mathrm{\β}}_k,T_s)j_{k-P}\left(t\right)+n\left(t\right)---\left(11\right)$

In formula (11), s (t) and j (t) represents the satellite-signal and undesired signal that receive, a (α respectively _k, β _k, T _s) be the steering vector of a kth echo signal (satellite-signal or undesired signal).T _stime domain delay line interval, α _kand β _kbe expressed as position angle and the angle of pitch that a kth echo signal arrives array antenna.N (t) represents white Gaussian noise, and its power spectrum density is expressed as N ₀/ 2.

In formula (11), a (α _k, β _k, T _s) represent space-time two-dimensional target vector, i.e. time vector a _s(T _s) and direction in space vector a _s(α _k, β _k) Crow Neck long-pending (KroneckerProduct), and can be expressed as:

$a({\mathrm{\α}}_k,{\mathrm{\β}}_k,T_s)=a({\mathrm{\α}}_k,{\mathrm{\β}}_k)\⊗a_s\left(T_s\right)---\left(12\right)$

The time delay unit number of each radio-frequency channel is made to be N, then

$a_s({\mathrm{\α}}_k,{\mathrm{\β}}_k)=\left[\begin{array}{cccc}1& e^{-j2\mathrm{\πf}{\mathrm{\τ}}_1}& \·\·\·& e^{-j2\mathrm{\πf}{\mathrm{\τ}}_{M-1}}\end{array}\right]---\left(13\right)$

$a_t=\left[\begin{array}{cccc}1& e^{-j2\mathrm{\πfT}}& \·\·\·& e^{-j2\mathrm{\πf}(N-1)T_s}\end{array}\right]---\left(14\right)$

In formula (13), f represents satellite frequency, τ _i(i=1,2 ... M-1) provided by formula (10).In formula (14), T _srepresent the time interval of delay cell, its value should be less than signal bandwidth.

Step 13: calculate Antenna array weight vector

Adopt multiple constraint minimum variance space-time adaptive process (MCMV-STAP) algorithm to retrain satellites in view signal, then make array output power minimum, thus suppress undesired signal while protection satellite-signal.This algorithm needs the best initial weights of trying to achieve array, is expressed as follows:

$\left\{\begin{array}{c}\mathrm{\ω}=\arg \min {\mathrm{\ω}}^H{R}_U\mathrm{\ω}\\ \mathrm{subjectto}:{\mathrm{\ω}}^{H}A=F^T\end{array}\right.---\left(15\right)$

In formula (15), ω represents the weighted vector of array, the conjugate transpose of H representing matrix, R _ube the array covariance matrix of input signal, can be expressed as

R _U=E{UU ^H}(16)

In formula (15), A represents the constraint matrix of satellite-signal, and F represents the constraint vector corresponding with A, A and F is expressed as

A=[a(α ₁,β ₁,T _s),a(α ₂,β ₂,T _s)…a(α _P,β _P,T _s)](17)

F ^t=[1,1 ... 1] _{1 × p}(18) adopt method of Lagrange multipliers, array weight vector ω can be expressed as

$\mathrm{\ω}=R_U^{-1}A{\left(A^H{R}_U^{-1}A\right)}^{-1}f---\left(19\right)$

Acquisition array weight after, can obtain array output express be for

y(t)=ω ^HU(t)。

Traditional space-time adaptive processing method, usually adopts Blind adaptive beamforming algorithm, does not namely need to know that satellite arrives position angle and the angle of pitch of array antenna, just simple a certain bay is carried out weights constraint as with reference to array element.But, owing to not retraining in satellite-signal direction, wave beam can not be formed in satellite-signal direction.Therefore, such blind beamforming algorithm also will weaken the satellite-signal expected while suppressing undesired signal, and satellite-signal power will be weakened, and can not obtain the maximize SINR of signal.

The carrier positions that the present invention mainly utilizes GPS/SINS integrated navigation system to provide and attitude, the satellite position that GPS ephemeris provides, thus obtain the steering vector between satellite to carrier, namely satellite arrives deflection and the angle of pitch of carrier.Then, this steering vector as the prior imformation of multiple constraint minimum variance space-time adaptive process (MCMV-STAP) algorithm, in spatial domain, time domain suppress simultaneously broadband interference and arrowband interference.

Steering vector method between calculating satellite to carrier of the present invention, it is characterized in that GPS/SINS calculate carrier position be arranged in earth coordinates, attitude is arranged in sky, northeast coordinate system, the satellite position utilizing satellite ephemeris to resolve is in ECEF coordinate system, and the steering vector needed for Wave beam forming should provide under carrier coordinate system.Tried to achieve prior imformation (position of satellite, the position of carrier and attitude) after changes in coordinates, carrier coordinate system is obtained deflection and the angle of pitch that satellite arrives carrier by the present invention.

The mode that the present invention adopts array antenna to combine with GPS/SINS is to improve the antijamming capability of navigation neceiver.In array antenna, adopt the space-time adaptive beamforming algorithm (MCMV-STAP) of multiple constraint minimum variance principle, this algorithm can form wave beam in many satellites in view directions, forms zero and falls into, thus suppress undesired signal while strengthening satellite-signal at interference radiating way.In order to obtain the peak power of satellite-signal, this algorithm needs to know that satellite-signal arrives position angle and the angle of pitch of carrier.Therefore, the present invention adopts circular array antenna structure, the Wave beam forming that introducing GPS/SINS integrated navigation is array antenna provides position and the attitude of carrier, satellite ephemeris is adopted to provide the position of satellite, thus provide prior imformation for Wave beam forming, namely satellites in view arrives position angle and the angle of pitch of carrier.

Accompanying drawing explanation

Fig. 1 is the implementing procedure figure of the associating anti-interference method based on array antenna and GPS/SINS;

Fig. 2 is the schematic diagram of array antenna, sky, northeast coordinate system and carrier coordinate system;

Fig. 3 is the structural drawing of multiple constraint minimum variance space-time adaptive processing (MCMV-STAP).

Embodiment

Below in conjunction with accompanying drawing, the present invention is further elaborated.

Step 1: initialization carrier positions, speed and attitude information

As shown in Figure 1, position and the attitude of first initialization carrier is needed.In earth coordinates, the coordinate of setting carrier initial time: latitude L=45 °, longitude λ=126 ° and height h=200m.The attitude angle of initialization carrier, comprises pitching angle theta=0 °, roll angle γ=0 ° and ψ=45 °, position angle.Then set the flight path of carrier, static, the track such as rectilinear motion or circular motion can be set as.The east orientation speed V of carrier flight path setting carrier as requested _e, north orientation speed V _nwith sky to speed V _u.Thus obtain the output of desirable gyroscope and accelerometer, f _e, f _nand f _ube expressed as accelerometer in east orientation, north orientation and sky to output.

Step 2: the parameter information of initialization Navigation Filter

In GPS/SINS integrated navigation wave filter, adopt feedback compensation mode.The quantity of state of Navigation Filter is attitude error, velocity error, site error, gyro error and accelerometer error.At initial time, carrier angle of pitch error φ _e=0.1 °, roll angle error φ _n=0.1 ° and azimuth angle error φ _u=1 °.The east orientation velocity error δ V of carrier _e=0.2m/s, north orientation velocity error δ V _n=0.2m/s and sky are to velocity error δ V _u=0.5m/s.Latitude error δ L=0 ° of carrier, longitude error δ λ=0 ° and height error δ h=20m.Gyro drift ε _b=0.1 °/h, white noise ε _g=0.05 °/h, accelerometer constant error white noise w _a=5 × 10 ^-4g, g are acceleration of gravity.Gyro error ε in east orientation, north orientation and sky to drift be expressed as ε _e, ε _nand ε _u.Accelerometer in east orientation, north orientation and sky to output error can be expressed as respectively with

Step 3: the error rate of computing gyroscope and accelerometer

Gyrostatic constant value drift can be described as with arbitrary constant

${\stackrel{\·}{\mathrm{\ϵ}}}_b=0---\left(1\right)$

Gyrostatic noise ε _gcan describe with Dirac function.Therefore, gyro error rate of change can be expressed as

$\stackrel{\·}{\mathrm{\ϵ}}={\stackrel{\·}{\mathrm{\ϵ}}}_b+{\mathrm{\ϵ}}_g---\left(2\right)$

For accelerometer error rate of change can be thought of as first-order Markov process, error model is taken as

$\stackrel{\·}{\▿}=-\frac{1}{T_a}\▿+w_a---\left(3\right)$

Wherein, T _arepresent correlation time, w _afor white-noise process.

Step 4: according to the parameter in step 1-3, calculates attitude of carrier angle error rate of change, velocity error rate of change and site error rate of change

1, attitude of carrier angle error rate of change is calculated

with be expressed as the rate of change of carrier angle of pitch error, roll angle error and azimuth angle error, attitude of carrier angle error rate of change can be expressed as

$\begin{array}{c}{\stackrel{\·}{\mathrm{\φ}}}_E=-\frac{{\mathrm{\δV}}_N}{R+h}+({\mathrm{\ω}}_{\mathrm{ie}}\sin L+\frac{V_E}{R+h}\tan L){\mathrm{\φ}}_N-({\mathrm{\ω}}_{\mathrm{ie}}\cos L+\frac{V_E}{R+h}\tan L){\mathrm{\φ}}_U+{\mathrm{\ϵ}}_E\\ {\stackrel{\·}{\mathrm{\φ}}}_N=\frac{{\mathrm{\δV}}_E}{R+h}-{\mathrm{\ω}}_{\mathrm{ie}}\sin \mathrm{L\δL}-({\mathrm{\ω}}_{\mathrm{ie}}\sin L+\frac{V_E}{R+h}\tan L){\mathrm{\φ}}_E-\frac{V_E}{R+h}{\mathrm{\φ}}_U+{\mathrm{\ϵ}}_N\\ {\stackrel{\·}{\mathrm{\φ}}}_U=\frac{{\mathrm{\δV}}_E}{R+h}\tan L+({\mathrm{\ω}}_{\mathrm{ie}}\cos L+\frac{V_E}{R+h}{\sec }^{2}L)\mathrm{\δL}+({\mathrm{\ω}}_{\mathrm{ie}}\cos L+\frac{V_E}{R+h}){\mathrm{\φ}}_E+\frac{V_N}{R+h}{\mathrm{\φ}}_N+{\mathrm{\ϵ}}_U\end{array}---\left(4\right)$

In formula (4), R represents earth radius, ω _ierepresent earth rotation angular speed.

2, the rate of change of bearer rate error is calculated

with be expressed as the east orientation velocity error of carrier, north orientation velocity error and the sky rate of change to velocity error, bearer rate error rate can be expressed as

$\begin{array}{c}{\stackrel{\·}{\mathrm{\δV}}}_E=(\frac{V_N}{R+h}\tan L-\frac{V_U}{R+h}){\mathrm{\δV}}_E+({2\mathrm{\ω}}_{\mathrm{ie}}\sin L+\frac{V_E}{R+h}\tan L){\mathrm{\δV}}_N\\ -(-{\mathrm{\ω}}_{\mathrm{ie}}\cos L+\frac{V_E}{R+h}){\mathrm{\δV}}_U+({2\mathrm{\ω}}_{\mathrm{ie}}{V}_N\cos L+\frac{V_E{V}_N}{R+h}{\sec }^{2}L+2{\mathrm{\ω}}_{\mathrm{ie}}{V}_U\sin L)\mathrm{\δL}\\ +{\▿}_E+f_N{\mathrm{\φ}}_U-f_U{\mathrm{\φ}}_N\\ {\stackrel{\·}{\mathrm{\δV}}}_N=-2({\mathrm{\ω}}_{\mathrm{ie}}\sin L+\frac{V_E}{R+h}\tan L){\mathrm{\δV}}_E-\frac{V_U}{R+h}{\mathrm{\δV}}_N-\frac{V_N}{R+h}{\mathrm{\δV}}_U\\ -(2{\mathrm{\ω}}_{\mathrm{ie}}{V}_E\cos L+\frac{V_E^2{\sec }^{2}L}{R+h})\mathrm{\δL}+{\▿}_N-f_E{\mathrm{\φ}}_U+f_U{\mathrm{\φ}}_E\\ {\stackrel{\·}{\mathrm{\δV}}}_U=(2{\mathrm{\ω}}_{\mathrm{ie}}\cos L+\frac{V_E}{R+h}){\mathrm{\δV}}_E+\frac{{2V}_N}{R+h}{\mathrm{\δV}}_N-2{\mathrm{\ω}}_{\mathrm{ie}}{V}_E\sin \mathrm{L\δL}+{\▿}_U+f_E{\mathrm{\φ}}_N-f_N{\mathrm{\φ}}_E\end{array}---\left(5\right)$

3, the rate of change of carrier positions error is calculated

with be expressed as the rate of change of the latitude error of carrier, longitude error and height error.Carrier positions error rate can be expressed as

$\begin{array}{c}\mathrm{\δ}\stackrel{\·}{L}=\frac{{\mathrm{\δV}}_N}{R+h}\\ \mathrm{\δ}\stackrel{\·}{\mathrm{\λ}}=\frac{{\mathrm{\δV}}_E}{R+h}\\ \mathrm{\δ}\stackrel{\·}{h}=\mathrm{\δ}{V}_U\end{array},\sec L+\frac{V_E}{R+h}\sec \mathrm{LtgL\δL}---\left(6\right)$

Step 5: introduce GPS pseudo range measurement information, adopts feedback compensation mode to correct SINS output information, obtains current

Position and attitude accurately.

Code phase error measured by GPS track loop calculates the pseudorange between satellite to carrier, pseudorange upgrades the quantity of state of Navigation Filter as measurement information, thus dope site error (the δ L that current SINS measures, δ λ, δ h) and attitude error (δ θ, δ γ, δ ψ).Then the position (L, λ, h) adopting feedback compensation mode to export SINS and attitude (θ, γ, ψ) correct, and obtain current position accurately and attitude.

L=L+δL

λ=λ+δλ

h=h+δh

θ=θ+δθ(7)

γ=γ+δγ

ψ=ψ+δψ

Step 6: calculate the coordinate of carrier in ECEF coordinate system

The carrier positions calculated by step 5 is arranged in earth coordinates, i.e. (L, λ, h), is translated in ECEF coordinate system, its position coordinates can be expressed as

$X_p^e={[(R+h)\cos L\cos \mathrm{\λ},(R+h)\cos L\sin \mathrm{\λ},(R+h)\sin L]}^T---\left(8\right)$

Step 7: calculate carrier coordinate system (b system) is tied to ECEF coordinate system (e system) transition matrix to the transition matrix of navigational coordinate system (n system), navigation coordinate

The attitude of carrier calculated by step 5 can calculate the transition matrix of carrier coordinate system to navigational coordinate system for

$C_b^n=\left[\begin{array}{ccc}\cos \mathrm{\γ}\cos \mathrm{\ψ}-\sin \mathrm{\ψ}\sin \mathrm{\θ}\sin \mathrm{\γ}& -\sin \mathrm{\ψ}\cos \mathrm{\θ}& \sin \mathrm{\γ}\cos \mathrm{\ψ}+\cos \mathrm{\γ}\sin \mathrm{\θ}\sin \mathrm{\ψ}\\ \cos \mathrm{\γ}\sin \mathrm{\ψ}+\sin \mathrm{\γ}\sin \mathrm{\θ}\cos \mathrm{\ψ}& \cos \mathrm{\ψ}\cos \mathrm{\θ}& \sin \mathrm{\γ}\sin \mathrm{\ψ}-\cos \mathrm{\γ}\sin \mathrm{\θ}\cos \mathrm{\ψ}\\ -\sin \mathrm{\γ}\cos \mathrm{\θ}& \sin \mathrm{\θ}& \cos \mathrm{\θ\γ}\cos \end{array}\right]---\left(9\right)$

The carrier positions calculated by step 5, can calculate the transition matrix that navigation coordinate is tied to ECEF coordinate system for

$C_n^e=\left[\begin{array}{ccc}-\sin \mathrm{\λ}& -\sin L\cos \mathrm{\λ}& \cos L\cos \mathrm{\λ}\\ \cos \mathrm{\λ}& -\sin L\sin \mathrm{\λ}& \cos L\sin \mathrm{\λ}\\ 0& \cos L& \sin L\end{array}\right]---\left(10\right)$

Step 8: calculate the steering vector between carrier to satellite

Can in the hope of the transformed matrix of carrier coordinate system to ECEF coordinate system by (9), (10) in step 7

$C_b^e=C_n^e{C}_b^n---\left(11\right)$

Satellite ephemeris is utilized to calculate the position of satellite in ECEF coordinate system and in conjunction with formula (8) (11), can in the hope of the steering vector between carrier to satellite

${\stackrel{\→}{r}}^b={\left(C_b^e\right)}^T(X_s^e-X_p^e)---\left(12\right)$

Step 9: in carrier coordinate system, calculates position angle and the angle of pitch that satellite arrives antenna

If will coordinate definition in carrier coordinate system is then satellite arrives the azimuth angle alpha of antenna and angle of pitch β can be expressed as respectively

$\mathrm{\α}=\arctan (\hat{x},\hat{y})---\left(13\right)$

$\mathrm{\β}=\arctan (\hat{z},\sqrt{{\hat{x}}^2+{\hat{y}}^2})---\left(14\right)$

Step 10: design circular array antenna structure

In order to control beam position satellites in view direction at position angle and angle of pitch direction simultaneously, the present invention adopts circular array antenna structure, as shown in Figure 2.6 array element is uniformly distributed in round battle array circumferentially, and make radius of circle be r, then the interval circumferentially between adjacent array element is also r.Choosing of array element interval, and to meet Nyquist's theorem the same at time-domain sampling interval, and Space domain sampling interval d should be less than 1/2 of satellite carrier wavelength X.By satellite frequency f=1575.42 × 10 ⁶mHz, so array element distance is

$d\≤\frac{\mathrm{\λ}}{2}=\frac{1}{2}\·\frac{c}{f_{l1}}=\frac{1}{2}\·\frac{3\×{10}^8}{1575.42\×{10}^6}=\frac{1}{2}\·0.19=0.095m=9.5\mathrm{cm}---\left(15\right)$

In formula (15), c is the light velocity.

In order to make main beam width narrower, secondary lobe is lower, and resolution is high, and what need to be tried one's best in array element interval is large, so get radius of circle r=d=9.5cm, the diameter of whole array antenna is about 19cm.

With antenna array place plane for xoy plane, take true origin as reference, x-axis array element is No. 1 array element, then the polar coordinates of 6 array elements are respectively: (r, 0), (r, π/3), (r, 2 π/3), (r, π), (r, 4 π/3), (r, 5 π/3).

Step 11: calculate the time delay between each array element of satellite-signal arrival antenna

According to (13) (14) in step 9, can by α and β representation unit vector

e(α,β)=(sinαcosβ,cosαcosβ,sinβ) ^T(16)

Therefore, satellite-signal arrives i-th bay and the mistiming τ i between first reference array element can be expressed as

τ _i=e ^t(x _i-x ₁)/ci=1,2 ... in M-1 (17) formula (17), M representative antennas array element number.

Step 12: set up array antenna received signals model

User generally can receive the satellite-signal of more than 4, therefore, needs to form satellite corresponding to multiple beam position.Suppose that array antenna received has arrived P satellite-signal, Q undesired signal, then the signal model that antenna receives can be expressed as

$U\left(t\right)=\underset{k=1}{\overset{P}{\mathrm{\Σ}}}a({\mathrm{\α}}_k,{\mathrm{\β}}_k,T_s)s_k\left(t\right)+\underset{k=P+1}{\overset{P+Q}{\mathrm{\Σ}}}a({\mathrm{\α}}_k,{\mathrm{\β}}_k,T_s)j_{k-P}\left(t\right)+n\left(t\right)---\left(18\right)$

In formula (18), s (t) and j (t) represents the satellite-signal and undesired signal that receive, a (α respectively _k, β _k, T _s) be the steering vector of a kth echo signal (satellite-signal or undesired signal).T _stime domain delay line interval, α _kand β _kbe expressed as position angle and the angle of pitch that a kth echo signal arrives array antenna.N (t) represents white Gaussian noise, and its power spectrum density is expressed as N ₀/ 2.

In formula (18), a (α _k, β _k, T _s) represent space-time two-dimensional target vector, i.e. time vector a _s(T _s) and direction in space vector a _s(α _k, β _k) Crow Neck long-pending (KroneckerProduct), and can be expressed as:

$a({\mathrm{\α}}_k,{\mathrm{\β}}_k,T_s)=a_s({\mathrm{\α}}_k,{\mathrm{\β}}_k)\⊗a_s\left(T_s\right)---\left(19\right)$

The time delay unit number of each radio-frequency channel is made to be N, then

$a_s({\mathrm{\α}}_k,{\mathrm{\β}}_k)=\left[\begin{array}{cccc}1& e^{-j2\mathrm{\πf}{\mathrm{\τ}}_1}& \·\·\·& e^{-j2\mathrm{\πf}{\mathrm{\τ}}_{M-1}}\end{array}\right]---\left(20\right)$

$a_t=\left[\begin{array}{cccc}1& e^{-j2\mathrm{\πfT}}& \·\·\·& e^{-j2\mathrm{\πf}(N-1)T_s}\end{array}\right]---\left(21\right)$

In formula (20), f represents satellite frequency, τ _i(i=1,2 ... M-1) provided by formula (17).In formula (21), T _srepresent the time interval of delay cell, its value should be less than signal bandwidth.

Step 13: calculate Antenna array weight vector

Accompanying drawing 3 is the structural drawing of multiple constraint minimum variance space-time adaptive processing (MCMV-STAP).From the passage of each array element, time delay at different levels constitutes FIR filter, can remove interference in time domain; From identical time delay node, different array element constitutes the auto adapted filtering in spatial domain, can resolve spatial interference source and then form spatial domain zero and fall into and suppress interference from spatial domain.And the process in spatial domain also can utilize the feedback information after Time Domain Processing further, therefore space time processing also has the ability simultaneously rejecting interference in space-time two-dimensional territory.

Adopt MCMV-STAP algorithm to retrain satellites in view signal, then make array output power minimum, thus suppress undesired signal while protection satellite-signal.This algorithm needs the best initial weights of trying to achieve array, is expressed as follows:

$\left\{\begin{array}{c}\mathrm{\ω}=\arg \min {\mathrm{\ω}}^H{R}_U\mathrm{\ω}\\ \mathrm{subjectto}:{\mathrm{\ω}}^{H}A=F^T\end{array}\right.---\left(22\right)$

In formula (22), ω represents the weighted vector of array, the conjugate transpose of H representing matrix, R _ube the array covariance matrix of input signal, can be expressed as

R _U=E{UU ^H}(23)

In formula (23), A represents the constraint matrix of satellite-signal, and F represents the constraint vector corresponding with A, A and F is expressed as

A=[a(α ₁,β ₁,T _s),a(α ₂,β ₂,T _s)…a(α _P,β _P,T _s)](24)

F ^t=[1,1 ... 1] _{1 × p}(25) adopt method of Lagrange multipliers, array weight vector ω can be expressed as

$\mathrm{\ω}=R_U^{-1}A{\left(A^H{R}_U^{-1}A\right)}^{-1}f---\left(26\right)$

Acquisition array weight after, can obtain array output express be for

y(t)=ω ^HU(t)。(27)。

Claims (1)

1., based on an associating anti-interference method of array antenna and GPS/SINS, after the position of initialization carrier and attitude, set up GPS/SINS integrated navigation state equation and measurement equation; GPS/SINS integrated navigation provides position and the attitude of carrier in real time, calculates the position of present satellites according to satellite ephemeris information simultaneously, thus obtains the steering vector between satellite to carrier, and namely satellite arrives deflection and the angle of pitch of carrier; Then, described steering vector as the prior imformation of multiple constraint minimum variance space-time adaptive Processing Algorithm, in spatial domain, time domain suppress simultaneously broadband interference and arrowband interference; It is characterized in that specifically comprising the steps:

Step 1: initialization carrier positions, speed and attitude information;

In earth coordinates, the coordinate of setting carrier initial time: latitude L, longitude λ and height h; The speed of initialization carrier northeastward in sky coordinate system: east orientation speed V _e, north orientation speed V _nwith sky to speed V _u; The attitude angle of initialization carrier, comprises pitching angle theta, roll angle γ and position angle ψ; Then the flight path setting carrier is the one in static, rectilinear motion or circular motion track; Thus obtain the output of gyroscope and accelerometer, f _e, f _nand f _ube expressed as accelerometer in east orientation, north orientation and sky to output;

Step 2: the parameter information of initialization Navigation Filter;

In GPS/SINS integrated navigation wave filter, adopt feedback compensation mode; The quantity of state of Navigation Filter is attitude error, velocity error, site error, gyro error and accelerometer error, φ _e, φ _nand φ _ube expressed as carrier angle of pitch error, roll angle error and azimuth angle error, δ V _e, δ V _nwith δ V _ube expressed as the east orientation velocity error of carrier, north orientation velocity error and sky to velocity error, δ L, δ λ and δ h are expressed as the latitude error of carrier, longitude error and height error, ε _e, ε _nand ε _ube expressed as gyroscope in east orientation, north orientation and sky to drift, ▽ _e, ▽ _nand ▽ _ube expressed as accelerometer in east orientation, north orientation and sky to output error;

Step 3: the error rate of computing gyroscope and accelerometer;

Step 4: according to the parameter in step 1-3, calculates attitude of carrier angle error rate of change, velocity error rate of change and site error rate of change;

Step 5: introduce GPS pseudo range measurement information, adopts feedback compensation mode to correct SINS output information, obtains current position accurately and attitude;

Code phase error measured by GPS track loop calculates the pseudorange between satellite to carrier, pseudorange upgrades the quantity of state of Navigation Filter as measurement information, thus dope the error of current SINS institute state quantity measurement amount, then adopt feedback compensation mode to correct SINS output information, obtain current position accurately and attitude;

Step 6: calculate the coordinate of carrier in ECEF coordinate system;

The carrier positions calculated by step 5 is arranged in earth coordinates, i.e. (L, λ, h), is translated in ECEF coordinate system, its position coordinates be expressed as

$X_p^e={\[(R+h)\cos Lcos\mathrm{\λ},(R+h)\cos Lsin\mathrm{\λ},(R+h)\sin L\]}^T$

Wherein, R is earth radius;

Step 7: calculate the transition matrix that carrier coordinate system and b are tied to the transition matrix of navigational coordinate system and n system, navigation coordinate is tied to ECEF coordinate system and e system;

The attitude of carrier angle calculated by step 5, calculates the transition matrix of carrier coordinate system to navigational coordinate system for

$C_b^n=\left[\begin{array}{ccc}\cos \mathrm{\γ}\cos \mathrm{\ψ}-\sin \mathrm{\ψ}\sin \mathrm{\θ}\sin \mathrm{\γ}& -\sin \mathrm{\ψ}\cos \mathrm{\θ}& \sin \mathrm{\γ}\cos \mathrm{\ψ}+\cos \mathrm{\γ}\sin \mathrm{\θ}\sin \mathrm{\ψ}\\ \cos \mathrm{\γ}\sin \mathrm{\ψ}+\sin \mathrm{\γ}\sin \mathrm{\θ}\cos \mathrm{\ψ}& \cos \mathrm{\ψ}\cos \mathrm{\θ}& \sin \mathrm{\γ}\sin \mathrm{\ψ}-\cos \mathrm{\γ}\sin \mathrm{\θ}\cos \mathrm{\ψ}\\ -\sin \mathrm{\γ}\cos \mathrm{\θ}& \sin \mathrm{\θ}& \cos \mathrm{\θ}\cos \mathrm{\γ}\end{array}\right]$

The carrier positions calculated by step 5, calculates the transition matrix that navigation coordinate is tied to ECEF coordinate system for

$C_n^e=\left[\begin{array}{ccc}-\sin \mathrm{\λ}& -\sin L\cos \mathrm{\λ}& \cos L\cos \mathrm{\λ}\\ \cos \mathrm{\λ}& -\sin L\cos \mathrm{\λ}& \cos L\sin \mathrm{\λ}\\ 0& \cos L& \sin L\end{array}\right];$

Step 8: calculate the steering vector between carrier to satellite;

By the transition matrix in step 7 transition matrix try to achieve the transformed matrix of carrier coordinate system to ECEF coordinate system

$C_b^e=C_n^e{C}_b^n$

Satellite ephemeris is utilized to calculate the position of satellite in ECEF coordinate system and in conjunction with formula $X_p^e={\[(R+h)\cos Lcos\mathrm{\λ},(R+h)\cos L\sin \mathrm{\λ},(R+h)\sin L\]}^T$ With $C_b^e=C_n^e{C}_b^n,$ Try to achieve the steering vector between carrier to satellite

${\stackrel{\→}{r}}^b={\left(C_b^e\right)}^T(X_s^e-X_p^e)$

Step 9: in carrier coordinate system, calculates position angle and the angle of pitch that satellite arrives antenna;

If will coordinate definition in carrier coordinate system is then satellite arrives the azimuth angle alpha of antenna and angle of pitch β can be expressed as respectively

$\mathrm{\α}=arctan(\hat{x},\hat{y})$

$\mathrm{\β}=arctan(\hat{z},\sqrt{{\hat{x}}^2+{\hat{y}}^2})$

Step 10: design circular array antenna structure;

6 array element is uniformly distributed in round battle array circumferentially, radius of circle is made to be r, interval then circumferentially between adjacent array element is also r, choosing of array element interval, and it is the same that time-domain sampling interval meets Nyquist's theorem, Space domain sampling interval d should be less than 1/2 of satellite carrier wavelength X, by satellite frequency f=1575.42 × 10 ⁶mHz, so array element distance is

$d\≤\frac{\mathrm{\λ}}{2}=\frac{1}{2}\·\frac{c}{f_{l1}}=\frac{1}{2}\·\frac{3\×{10}^8}{1575.42\×{10}^6}=\frac{1}{2}\·0.19=0.095m=9.5cm$

Wherein, c is the light velocity;

Step 11: calculate the time delay between each array element of satellite-signal arrival antenna;

According in step 9 $\mathrm{\α}=arctan(\hat{x},\hat{y}),$ $\mathrm{\β}=arctan(\hat{z},\sqrt{{\hat{x}}^2+{\hat{y}}^2}),$ By α and β representation unit vector

e(α,β)＝(sinαcosβ,cosαcosβ,sinβ) ^T

Therefore, satellite-signal arrives i-th bay and to the mistiming τ between first reference array element _ican be expressed as

τ _i＝e ^T·(x _i-x ₁)/ci＝1,2,…M-1

Wherein, M representative antennas array element number;

Step 12: set up array antenna received signals model

If array antenna received has arrived P satellite-signal, Q undesired signal, then the signal model that antenna receives has been expressed as

$U\left(t\right)=\underset{k=1}{\overset{P}{\Σ}}a({\mathrm{\α}}_k,{\mathrm{\β}}_k,T_s)s_k\left(t\right)+\underset{k=P+1}{\overset{P+Q}{\Σ}}a({\mathrm{\α}}_k,{\mathrm{\β}}_k,T_s)j_{k-P}\left(t\right)+n\left(t\right)$

Wherein, s (t) and j (t) represents the satellite-signal and undesired signal that receive, a (α respectively _k, β _k, T _s) be the steering vector of a kth echo signal, T _stime domain delay line interval, α _kand β _kbe expressed as position angle and the angle of pitch that kth echo signal arrives array antenna, n (t) represents white Gaussian noise, its power spectrum density is expressed as N ₀/ 2;

A (α _k, β _k, T _s) represent space-time two-dimensional target vector, i.e. time vector a _s(T _s) and direction in space vector a _s(α _k, β _k) Crow Neck amass, and to be expressed as:

$a({\mathrm{\α}}_k,{\mathrm{\β}}_k,T_s)=a_s({\mathrm{\α}}_k,{\mathrm{\β}}_k)\⊗a_s\left(T_s\right)$

The time delay unit number of each radio-frequency channel is made to be N, then

$a_s({\mathrm{\α}}_k,{\mathrm{\β}}_k)=\left[\begin{array}{cccc}1& e^{-j2{\mathrm{\πf\τ}}_1}& \mathrm{...}& e^{-j2{\mathrm{\πf\τ}}_{M-1}}\end{array}\right]$

a _t＝[1e ^-j2πfT…e ^{-j2πf(N-1)Ts}]

F represents satellite frequency, T _srepresent that the time interval, its value of delay cell should be less than signal bandwidth;

Step 13: calculate Antenna array weight vector;

Adopt multiple constraint minimum variance space-time adaptive Processing Algorithm to retrain satellites in view signal, then make array output power minimum, being expressed as follows of the best initial weights of array:

$\left\{\begin{array}{c}\mathrm{\ω}=\arg {\mathrm{min\ω}}^H{R}_U\mathrm{\ω}\\ subjectto:{\mathrm{\ω}}^{H}A=F^T\end{array}\right.$

Wherein, ω represents the weighted vector of array, the conjugate transpose of H representing matrix, R _ube the array covariance matrix of input signal, be expressed as

R _U＝E{UU ^H}

A represents the constraint matrix of satellite-signal, and F represents the constraint vector corresponding with A, A and F is expressed as

A＝[a(α ₁,β ₁,T _s),a(α ₂,β ₂,T _s)…a(α _P,β _P,T _s)]

f ^T＝[1,1…1] _1×p

Adopt method of Lagrange multipliers, array weight vector ω is expressed as

$\mathrm{\ω}=R_U^{-1}A{\left(A^H{R}_U^{-1}A\right)}^{-1}f$

Acquisition array weight after, obtain array output express be for

y(t)＝ω ^HU(t)。