CN104155662B - The mutual disturbance restraining method of self adaptation based on GNSS correlation peak detector - Google Patents

Abstract

The mutual disturbance restraining method of self adaptation based on GNSS correlation peak detector, utilizes navigation neceiver to receive desired signal and interference signal simultaneously, obtains disturbing signal to affect relation to what desired signal captured.Then at desired signal and the Doppler frequency shift Δ f of interference signal_iFor setting up grid on the code delay between two adjacent integers times of KHz and Doppler frequency shift, set up based on χ for each grid²Distribution differentiates that the GNSS of hypothesis testing disturbs correlation peak detector mutually.If detecting cross-correlation interference, then scanned for by interference signal mutual to GNSS, construct corresponding interference space, self adaptation subspace projection algorithm is utilized to carry out AF panel, and utilize the mutual interference mitigation gains judgment criterion as error of output chip segment signal, determine whether to obtain desired signal by error judgment.If error is excessive does not obtains desired signal, the most again carries out subspace projection and carry out AF panel, until obtaining the desired signal that disclosure satisfy that requirement.

Description

The mutual disturbance restraining method of self adaptation based on GNSS correlation peak detector

Technical field

The invention belongs to the communications field, relate to the mutual interfering signal suppressing method of a kind of satellite navigation system.

Background technology

GLONASS (GNSS) is various countries for military or civilian purpose, and develop a set of Satellite is used to provide the system of position and time.From guided missile, opportunity of combat and warship to automobile, civil aircraft, individual People's computer or even handheld communication devices, almost can use Technique of Satellite Navigation and Positioning.

Along with development and the continuous expansion of application thereof of Satellite Navigation Technique, the most numerous countries are numerous and confused Step up to build independent satellite navigation system.The U.S. is implementing GPS updating program, Russia Stepping up to recover GLONASS all round properties, European Union is building GALILEO, and China steps up to build COMPASS system.

At present, international telegraph union (ITU) exists exclusively in several frequency ranges of Radio Satellite navigation configuration Have the gps signal of the U.S. and Muscovite GLONASS signal, add plan launch GALILEO and COMPASS signal, and GPS Modernization Signal, these frequency ranges have become suitable crowded, except The frequency range of GLONASS does not has outside obvious overlap, and the frequency range that other three systems have partly overlaps, and has Completely overlapped, therefore, satellite navigation signals interference each other is inevitable.In general, same frequency In section, GNSS disturbs mutually and is divided into two kinds, from the signal phase of different satellites in the most same satellite navigation system The mutually reception of interfering with an opponent, such as interference etc. between different satellite navigation signals in GPS system, be referred to as from Interference or system internal interference；And the interference between two different satellite navigation systems, such as GPS system And the signal disturbing etc. between GALILEO system, referred to as inter-system interference.Additionally, utilizing satellite navigation During the reflected signal of system carries out remote sensing application, the direct signal of same satellite is for reflected signal Interference also referred to as self-interference.

In the application of navigation signal, in user terminal is in some special environments, such as in parking lot, Between architecture indoor and housing-group or in valley, countryside, forest and outlying district, different satellite-signals arrive and receive The dynamic range of the signal intensity of machine is very big, causes the signal power of different passage to have a long way to go.Therefore by strong The self-interference impact that weak signal cross-correlation causes, the weak signal under strong signal disturbing is difficult to capture, often sinks into by mistake Catch to lose the most completely and catch.Therefore, the detection of self-interference in same system and the research of rejecting are had very much Necessary.

According to investigation, for the detection of GNSS correlation peak, pertinent literature both domestic and external is considerably less, and skill Art itself is the most immature, does not has implementation method conveniently.Suppression after suffering strong cross-correlation interference Method, general employing counteracting serial interference method.The method detect interference signal waveform time, one time one Individual gradually interference is removed from receiving signal, concrete, can according to received signal power successively decrease suitable Sequence is demodulated, and the signal with strongest received power is the most demodulated, and detection demodulated at signal Afterwards, it is deducted from reception signal.But, at present in most GNSS receiver, receive The quantization digit of signal is only 2, and the method precision directly deducting interference from reception signal can be the lowest, it is difficult to Realize eliminating the purpose of high reject signal.

Summary of the invention

Present invention solves the technical problem that and be: overcome the deficiencies in the prior art, it is provided that be a kind of based on GNSS The mutual disturbance restraining method of self adaptation of correlation peak detector, not only solves satellite navigation system signals and does mutually Acquisition of signal problem when disturbing, additionally it is possible to the signal of interference mutually detected effectively is suppressed and rejects.

The technical solution of the present invention is: self adaptation based on GNSS correlation peak detector is disturbed mutually and pressed down Method processed, comprises the steps:

(1) utilize navigation neceiver to receive aeronautical satellite A and the navigation signal of aeronautical satellite B simultaneously, its The signal of middle aeronautical satellite A is desired signal, and the signal of aeronautical satellite B is interference signal；

(2) between Doppler frequency shift is KHz two adjacent integers times of desired signal and interference signal Code delay and Doppler frequency on set up lattice respectively, and Doppler frequency lattice and a code delay are divided Lattice merge one grid of formation, set up based on χ at each grid²Distribution differentiates hypothesis testing GNSS disturbs correlation peak detector mutually, utilizes GNSS mutual Interference Peaks value detector to carry out cross-correlation interference Detection, without detecting cross-correlation interference, then directly receive desired signal and terminate；If detection To cross-correlation interference, then proceed to step (3)；

Described GNSS disturbs the form of correlation peak detector to be mutually:

Wherein, X_i() is GNSS cross-correlation chi square function,

$X_i(\mathrm{\τ},\mathrm{\Δ}{f}_i)=\frac{1}{N_I}\underset{r=0}{\overset{N_I-1}{\mathrm{\Σ}}}{\left|{\∫}_{r{\mathrm{ZN}}_C}^{(r+1)Z{N}_C}x\left(t\right)c_i^{*}\left(R\left(\mathrm{\Δ}{f}_i\right)(t-\mathrm{\τ})\right)e^{-j2\mathrm{\π\Δ}{f}_{i}t}\mathrm{dt}\right|}^2$

N_cFor the number of coherent integration chip, N_IFor the number of coherent integration data block, ()^*Represent conjugate complex number, Ζ For chip-spaced length, r is integral number of times, and τ is code delay, and x (t) is desired signal, c_iT () is interference letter The chip value of number i, R (Δ f_i) it is that Doppler frequency shift is Δ f under preferable tracking condition_iTime correlation function value, m₀ For doppler cells, at m₀In unit, make X_i(n,m₀) n value that maximum correlation is corresponding is n₀, finding After maximum related value, find second largest correlation with sampling interval Δ, make X_i(n,m₀) second largest being correlated with The n value of value correspondence isFor the standard deviation estimate of correlation, TNF is detection thresholding,

N is the maximum of code phase dimension；

(3) utilize the weighter factor g of self adaptation subspace projection, the data signal received carried out cutting, It is divided into L the fragment signal that length is identical, for each section of fragment signal, does mutually according to what acquisition and tracking arrived Disturb the data bit information of signal, code phase, carrier frequency value and initial phase value, build and disturb letter mutually Number projection matrix, and utilize subspace projection algorithm, from data signal, remove interference signal, obtain each The residual signal of fragment signal enters step (4)；

WhereinX is the total length of the data signal received, and L is the fragment signal number being divided into Amount,Integer part, the value of K be numeral 1～20 in any integer；

(4) calculate the gain of the residual signal of each fragment signal respectively, and compare with gain threshold, right Return in step (3) after the gain fragment signal section of segmentation again beyond gain threshold and utilize subspace Projection algorithm removes interference signal, until the gain of the residual signal of new segmentation section is less than gain threshold；So After by original for the first time without departing from fragment signal and the section of segmentation and again segmenting again of gain threshold After Duan, each segmentation segment signal of the gain respectively less than gain threshold of each segmentation section is combined, it is thus achieved that desired signal.

The size of described Doppler frequency lattice be Doppler frequency shift be KHz two adjacent integers times it Between 2/3 times.1/2 chip that size is code phase of described code delay lattice.

Present invention advantage compared with prior art is:

(1) the inventive method for GNSS signal carrier phase and Doppler frequency shift the shadow to coherent signal The principle rung, it is proposed that based on χ²The cross-correlation interference peak detector of Testing Statistical Hypotheses carries out interference and visits Survey, and on this basis, construct corresponding interference space according to the interference signal parameters obtained, utilize son Space projection carries out AF panel, and the method had both reduced the error impact on jamming performance, made again amount of calculation Decline to a great extent, be suitable for changeable environment and many data samples situation；

(2) the inventive method is analyzed as a example by relevant DLL, draws at Doppler frequency shift Δ f_iFor KHz Integral multiple time, error that cross-correlation can produce is maximum, at Doppler frequency shift Δ f_iNon-integral multiple for KHz Time, the less conclusion of error that cross-correlation produces.Conclusion accordingly, when detecting cross-correlation interference, at KHz Two adjacent integers times between set up grid and test, there is typicality；And use χ²Distribution differentiates standard Then so that the detection of GNSS cross-correlation peak value is totally independent of interference noise so that the cross-correlation interference of foundation Peak detector has more versatility；

(3) the inventive method is in self adaptation subspace projection method, according to follow the tracks of strong signal amplitude, Phase place constructs the projection subspace of relative fragment with Doppler frequency shift, accordingly even when strong signal capture error code, right For subspace projection algorithm, occur that the data length of mistake only limits the fragment signal in an inferior subspace computing In, do not interfere with subspace projection computing next time, i.e. have no effect on the AF panel of whole data length.

Accompanying drawing explanation

Fig. 1 is the flow chart of the inventive method；

Fig. 2 is the schematic diagram of subspace projection.

Detailed description of the invention

As it is shown in figure 1, be the FB(flow block) of the inventive method, including following step:

Step 1: utilize navigation neceiver simultaneously to receive the navigation signal of aeronautical satellite A and other aeronautical satellite (signal assuming wherein aeronautical satellite A is desired signal, aeronautical satellite B or the signal of aeronautical satellite C For interference signal), obtain between the desired signal of aeronautical satellite A and the interference signal of other aeronautical satellite Cross-correlation interference, thus obtains disturbing signal to affect relation to what desired signal captured.

GNSS satellite navigation system is a kind of information system based on CDMA, the spread spectrum used due to it Code is not completely orthogonal, so its performance is faced with the impact of cross-correlation interference.GNSS receiver institute The cross-correlation interference faced is in the starting stage, this is because GNSS satellite distance terrestrial user geometry away from From all more than 20,000 kms, so the satellite-signal power difference that GNSS receiver receives is not very big, Cross-correlation interference is inconspicuous.But it is as the satellite navigation extreme environment such as in city and in building In application, when utilizing the signal that navigation neceiver receives different navigation satellite, due to the screening of building simultaneously Gear, between the power of the different satellite-signals that GNSS receiver receives, difference may be very big, desired signal Signal power may become the signal power much smaller than interference signal, now high reject signal and weak desired signal Between cross-correlation interference can directly influence the capture of desired signal.

The concept of GNSS satellite modulation and demodulation is to use different PRN code with each satellite, but has Based on having identical spreading rate and carrier wave, the technology of this modulation /demodulation is referred to as CDMA (CDMA).CDMA technology requires required the launched PRN of aeronautical satellite A of GNSS receiver synthesis The reproduction code of code, and the phase place of this reproduction code mobile makes the specific PRN of itself and every tracked satellite There is relevant (the such as signal of aeronautical satellite B or aeronautical satellite C is interference signal) in code.For whole code Any phase place or the combination of Doppler frequency shift in cycle, the PRN code of every satellite use is defended with any other The PRN code of star must have the cross-correlation of minimum.

The cross-correlation impact on GNSS signal for convenience of description, illustrates as a example by relevant DLL, Relevant DLL phase discriminator is:

$D\left[k\right]=({\tilde{Z}}_E{\left[k\right]}^2-{\tilde{Z}}_L{\left[k\right]}^2)/2---\left(1\right)$

Wherein: D [k] is the output error of k-th phase discriminator；

WithIt is respectively the correlation function of the lead and lag of k-th passage, is respectively as follows:

$\begin{array}{c}{\tilde{Z}}_E\left[k\right]={\tilde{S}}_E+\widetilde{\mathrm{\η}}\\ =\sqrt{2P_k}\exp \left(j{\mathrm{\Δ\θ}}_1\right)R_k(-\frac{{\mathrm{dT}}_c}{2})+\sqrt{2P_i}\exp \left(j{\mathrm{\Δ\θ}}_2\right)R_k^{\′}(-\frac{{\mathrm{dT}}_c}{2})+\widetilde{\mathrm{\η}}\end{array}---\left(2\right)$

$\begin{array}{c}{\tilde{Z}}_L\left[k\right]={\tilde{S}}_L+\widetilde{\mathrm{\η}}\\ =\sqrt{2P_k}\exp \left(j{\mathrm{\Δ\θ}}_1\right)R_k(+\frac{{\mathrm{dT}}_c}{2})+\sqrt{2P_i}\exp \left(j{\mathrm{\Δ\θ}}_2\right)R_k^{\′}(+\frac{{\mathrm{dT}}_c}{2})+\widetilde{\mathrm{\η}}\end{array}---\left(3\right)$

$R_k\left(\mathrm{\τ}\right)=\frac{1}{T_c}{\∫}_0^{T_c}{d}_k\left(t\right)d_k(t-\widehat{\mathrm{\τ}})c_k(t-\mathrm{\τ})c_k(t-\widehat{\mathrm{\τ}})\exp \left(j2\mathrm{\π\Δ}{f}_{k}t\right)\mathrm{dt}---\left(4\right)$

$R_k^{\′}\left(\mathrm{\τ}\right)=\frac{1}{T_c}{\∫}_0^{T_c}{d}_i(t-{\mathrm{\τ}}_i)d_k(t-\widehat{\mathrm{\τ}})c_i(t-{\mathrm{\τ}}_i)c_i(t-\widehat{\mathrm{\τ}})\exp \left(j2\mathrm{\π\Δ}{f}_i\right)\mathrm{dt}---\left(5\right)$

Under preferable tracking condition, it is assumed that entirely without the impact of cross-correlation, the phase of receiver k-th passage Pass functional value is:

$R_k\left(\frac{{\mathrm{dT}}_c}{2}\right)=R_k(-\frac{2{\mathrm{dT}}_c}{2})=\sqrt{2P_k}{T}_c(1-\frac{d}{2})---\left(6\right)$

Wherein, T_cFor coherent integration time, d is the real-time spacing of GNSS receiver correlator.For advanced relevant Function,For lag correlation function,For error, P_kFor the signal power of k-th passage, Δ θ₁It is K The phase value of individual channel signal, P_iFor disturbing the power of signal i, Δ θ₂For disturbing the phase value of signal i, d_k(t) The correlator distance values of k-th channel signal, d_iT the correlator distance values of () interference signal i, τ is GNSS Receiver code postpones,For the conjugate of τ, τ_iThe code delay value of interference signal i, c_k(t) k-th passage Signal code chip value, c_iThe chip value of (t) interference signal i, Δ f_kThe signal(-) carrier frequency Doppler of k-th passage Deviation, Δ f_iCarrier frequency Doppler's deviation of interference signal i, exp () is the expression-form of index.

By formula (6) it is known that under preferable tracking condition, D [k] is output as 0.

The PRN code that aeronautical satellite A uses and other satellite any (aeronautical satellite B or aeronautical satellite C) PRN code must cross-correlation, equal for the combination of any phase place or Doppler frequency shift within the whole code cycle It is such.Therefore, after signal capture, it is assumed that the Frequency Estimation obtained is substantially accurate, i.e. Δ f_k=0, and Auto-correlation functionWithCompare cross-correlation functionWithBig many, therefore, The output of k-th phase discriminator can be reduced to:

$\begin{array}{c}D\left[k\right]=({\tilde{Z}}_E{\left[k\right]}^2-{\tilde{Z}}_L{\left[k\right]}^2)/2=\left({\tilde{Z}}_E\right[k]-{\tilde{Z}}_L[k\left]\right)\left({\tilde{Z}}_E\right[k]+{\tilde{Z}}_L[k\left]\right)/2\\ \≅\sqrt{2P_i}\×\sqrt{2P_k}\×\exp \left({\mathrm{j\Δ\θ}}_2\right)\×(R_k^{\′}(-\frac{{\mathrm{dT}}_c}{2})-R_k^{\′}\left(\frac{d{T}_c}{2}\right))\×R_k\left(\frac{{\mathrm{dT}}_c}{2}\right)\\ =\sqrt{2P_i}\×\sqrt{2P_k}\×T_c\×(1-\frac{d}{2})\×\exp \left(j{\mathrm{\Δ\θ}}_2\right)\\ \×\underset{(k-1)T}{\overset{\mathrm{kT}}{\∫}}{d}_k(t-\widehat{\mathrm{\τ}})\·d_i(t-{\mathrm{\τ}}_i)\·\exp \left(2\mathrm{\π\Δ}{f}_{i}t\right)c_i(t-{\mathrm{\τ}}_i)\·[c_k(t-\frac{d}{2})-c_k(t+\frac{d}{2})]\mathrm{dt}\end{array}---\left(7\right)$

The time adding up the repeatedly integral result of GNSS correlator is referred to as the time of integration, when integration Between T more than 20ms time, code Doppler frequency shift is the least, on the relative delay road of each integration interval code Footpath ρ (t) is a constant, and within the time of integration, the change of the Doppler frequency shift of carrier wave is also the least, because of The Doppler frequency shift Δ f of this carrier wave_iIt is regarded as a constant, because τ_iIt is worth less, d_kAnd d_iValue can be "+1 " simultaneously Or "-1 ", then d_k(t)·d_i(t-τ_i)=1.

Owing to, in integration interval, the Doppler frequency shift of carrier wave and code delay are considered as a constant, therefore, Formula (7) is represented by:

$D\left[k\right]=\sqrt{2P_i}\×\sqrt{2P_k}\×T_c\×(1-\frac{d}{2})\×\exp \left(j{\mathrm{\Δ\θ}}_2\right)\×\underset{(k-1)T}{\overset{\mathrm{kT}}{\∫}}P\left(t\right)\exp \left(2\mathrm{\π\Δ}{f}_{i}t\right)\mathrm{dt}---\left(8\right)$

Wherein:

$P\left(t\right)=d_k\left(t\right)\·d_i(t-{\mathrm{\τ}}_i)\·c_i(t-{\mathrm{\τ}}_i)\·\left[c_k(t-\frac{d}{2})-c_k(t+\frac{d}{2})\right]---\left(9\right)$

When navigation data phase alignment, d_k(t)·d_i(t-τ_i)=1, P (t) may be expressed as:

$P\left(t\right)=c_i(t-{\mathrm{\τ}}_i)\·\left[c_k(t-\frac{d}{2})-c_k(t+\frac{d}{2})\right]---\left(10\right)$

Δf_iCan be broken down into two parts:

$\mathrm{\Δ}{f}_i=\frac{m_1}{T}\±m_2$

Wherein m₁For integer,

$D\left[k\right]=\sqrt{2P_i}\×\sqrt{2P_k}\×T_c\×(1-\frac{d}{2})\×\exp \left(j{\mathrm{\Δ\θ}}_2\right)\×\underset{(k-1)T}{\overset{\mathrm{kT}}{\∫}}P\left(t\right)\exp \left(2\mathrm{\π\Δ}{f}_{i}t\right)\mathrm{dt}---\left(8\right)$

Analyze the periodicity of D [k],

$\begin{array}{c}D[k+\frac{1}{m_{2}T}]\\ =\sqrt{2P_i}\×\sqrt{2P_k}\×T_c\×(1-\frac{d}{2})\×\exp \left(j{\mathrm{\Δ\θ}}_2\right)\×\underset{(k-1+\frac{1}{m_{2}T})T}{\overset{(k+\frac{1}{m_{2}T})T}{\∫}}P\left(t\right)\exp \left[2\mathrm{\π}(\frac{m_1}{T}\±m_2)t\right]\mathrm{dt}\\ =\sqrt{2P_i}\×\sqrt{2P_k}\×T_c\×(1-\frac{d}{2})\×\exp \left(j{\mathrm{\Δ\θ}}_2\right)\×\underset{(k-1)T}{\overset{\mathrm{kT}}{\∫}}P\left(t\right)\exp \left[2\mathrm{\π}(\frac{m_1}{T}\±m_2)(t+\frac{1}{m_2})\right]\mathrm{dt}\\ =\sqrt{2P_i}\×\sqrt{2P_k}\×T_c\×(1-\frac{d}{2})\×\exp \left(j{\mathrm{\Δ\θ}}_2\right)\×\underset{(k-1)T}{\overset{\mathrm{kT}}{\∫}}P\left(t\right)\exp \left\{\right[2\mathrm{\π}(\frac{m_1}{T}\±m_2)t]+2\mathrm{\π}(\frac{m_1}{T}\±m_2)\frac{1}{m_2}\}\mathrm{dt}\end{array}---\left(11\right)$

IfIt it is integer

$\begin{array}{c}D[k+\frac{1}{m_{2}T}]\\ =\sqrt{2P_i}\×\sqrt{2P_k}\×T_c\×(1-\frac{d}{2})\×\exp \left(j{\mathrm{\Δ\θ}}_2\right)\underset{(k-1)T}{\overset{\mathrm{kT}}{\∫}}P\left(t\right)\exp \left[2\mathrm{\π}(\frac{m_1}{T}+m_2)t\right]\mathrm{dt}\\ =D\left[k\right]\end{array}---\left(12\right)$

Therefore, ifBe integer, then it is the cycle of D [k].

M is set₃Value, it is assumed thatUtilize m₃The value of T judges the periodicity of D [k], then above formula is

$\begin{array}{c}D[k+\frac{1}{m_{3}T}]\\ =\sqrt{2P_i}\×\sqrt{2P_k}\×T_c\×(1-\frac{d}{2})\×\exp \left(j{\mathrm{\Δ\θ}}_2\right)\\ \×\underset{(k-1)T}{\overset{\mathrm{kT}}{\∫}}P\left(t\right)\exp \left[2\mathrm{\π}(\frac{m_1}{T}\±m_2)(t+\frac{1}{m_3})\right]\mathrm{dt}\\ =\sqrt{2P_i}\×\sqrt{2P_k}\×T_c\×(1-\frac{d}{2})\×\exp \left(j{\mathrm{\Δ\θ}}_2\right)\\ \×\underset{(k-1)T}{\overset{\mathrm{kT}}{\∫}}P\left(t\right)\exp \left\{\right[2\mathrm{\π}(\frac{m_1}{T}\±m_2)t]+2\mathrm{\π}\left(\frac{m_1}{T{m}_3}\right)\±2\mathrm{\π}\frac{m_2}{m_3}\}\mathrm{dt}\end{array}---\left(13\right)$

For any m₃>m₂,Above formula can not be equal to D [k], and therefore D [k] is a periodicity Function, andIt it is the minimum period of D [k].

$\begin{array}{c}D\left[k\right]=\sqrt{2P_i}\×\sqrt{2P_k}\×T_c\×(1-\frac{d}{2})\×\exp \left(j{\mathrm{\Δ\θ}}_2\right)\\ \×\underset{(k-1)T}{\overset{\mathrm{kT}}{\∫}}{c}_i(t-{\mathrm{\τ}}_i)\·[c_k(t-\frac{d}{2})-c_k(t+\frac{d}{2})\·\exp \left(2\mathrm{\π\Δ}{f}_{i}t\right)\mathrm{dt}\\ =\sqrt{2P_i}\×\sqrt{2P_k}\×T_c\×(1-\frac{d}{2})\×\exp \left(j{\mathrm{\Δ\θ}}_2\right)\\ \×\{\underset{(k-1)T}{\overset{\mathrm{kT}}{\∫}}{c}_i(t-{\mathrm{\τ}}_i)\·c_k(t-\frac{d}{2})\·\exp \left(2\mathrm{\π\Δ}{f}_{i}t\right)\mathrm{dt}-\underset{(k-1)T}{\overset{\mathrm{kT}}{\∫}}{c}_i(t-{\mathrm{\τ}}_i)\·c_k(t+\frac{d}{2})\·\left(2\mathrm{\π\Δ}{f}_{i}t\right)\mathrm{dt}\end{array}---\left(14\right)$

When standard correlators spacing d=1 chip, the cross-correlation function of lead and lag is completely self-contained, The subtracting each other part and can not deduct any cross-correlation of above formula, but, for little correlator spacing d, advanced and Delayed cross-correlation function connects each other, and their common ground will be deducted by offseting.

For T=1ms, for any Δ f_i≠ 0, cross-correlation can be maintained at certain level, for T > 1ms, The integration of above formula shows, at high Δ f_iIn the case of, cross-correlation energy can be reduced.

$\begin{array}{c}D\left[k\right]=\sqrt{2P_i}\×\sqrt{2P_k}\×T_c\×(1-\frac{d}{2})\×\exp \left(j{\mathrm{\Δ\θ}}_2\right)\\ \×\underset{0}{\overset{0.002}{\∫}}{c}_i(t-{\mathrm{\τ}}_i)\·[c_k(t-\frac{d}{2})-c_k(t+\frac{d}{2})]\·\exp \left(1000\mathrm{\πt}\right)\mathrm{dt}\\ =\sqrt{2P_i}\×\sqrt{2P_i}\×T_c\×(1-\frac{d}{2})\×\exp \left(j{\mathrm{\Δ\θ}}_2\right)\\ \×\{\underset{0}{\overset{0.002}{\∫}}{c}_i(t-{\mathrm{\τ}}_i)\·c_k(t-\frac{d}{2})\·\exp \left(1000\mathrm{\πt}\right)\mathrm{dt}-\underset{0}{\overset{0.002}{\∫}}{c}_i(t-{\mathrm{\τ}}_i)\·c_k(t+\frac{d}{2})\·\exp \left(1000\mathrm{\πt}\right)\mathrm{dt}\\ =\sqrt{2P_i}\×\sqrt{2P_k}\×T_c\×(1-\frac{d}{2})\×\exp \left(j{\mathrm{\Δ\θ}}_2\right)\\ \×\left\{\begin{array}{c}[\underset{0}{\overset{0.001}{\∫}}{c}_i(t-{\mathrm{\τ}}_i)\·c_k(t-\frac{d}{2})\·\exp \left(1000\mathrm{\πt}\right)\mathrm{dt}+\underset{0.001}{\overset{0.002}{\∫}}{c}_i(t-{\mathrm{\τ}}_i)\·\exp \left(1000\mathrm{\πt}\right)\mathrm{dt}]\\ +[\underset{0}{\overset{0.001}{\∫}}{c}_i(t-{\mathrm{\τ}}_i)\·c_k(t+\frac{d}{2})\·\exp \left(1000\mathrm{\πt}\right)\mathrm{dt}+\underset{0.001}{\overset{0.002}{\∫}}{c}_i(t-{\mathrm{\τ}}_i)\·c_k(t+\frac{d}{2})\·\exp \left(1000\mathrm{\πt}\right)\mathrm{dt}]\end{array}\right\}\end{array}---\left(15\right)$

From formula (15), as Δ f_iDuring=500Hz and T=2ms, CA code and carrier phase are at first 1ms The most contrary with second 1ms, the cross-correlation error of first 1ms and the cross-correlation of second 1ms It is 0 that error is cancelled out each other, accordingly, because Integral Processing is zero in lead and lag passage cross-correlation error, so And, as Δ f_iDuring=1KHz, CA code and carrier phase repeat second 1ms completely at first 1ms, because of This cross-correlation error is sum of the two, and Integral Processing can not eliminate any error.

It follows that to any time of integration, at Doppler frequency shift Δ f_iDuring for the integral multiple of KHz, mutually Close the error that can produce maximum, at Doppler frequency shift Δ f_iDuring for KHz non-integral multiple, the mistake that cross-correlation produces Difference is less, therefore, when detecting cross-correlation interference, sets up grid between two adjacent integers times of KHz Test, and interference value thresholding must be set when setting up cross-correlation interference peak detector.

Step 2: the problem that affects desired signal captured for aeronautical satellite cross-correlation interference, set up based on_χ2 Distribution differentiates the cross-correlation peak value detector of hypothesis testing.

The cross-correlation interference detection of GNSS signal is the process of a search, and interference signal reappears with receiver The code of satellite is relevant with carrier wave, and wherein GNSS receiver code delay τ associates with reproduction code-phase, and many Pu Lewei Δ f_iBeing associated with reproduction carrier wave, each code phase search increment is the lattice (typical case of a code Under, code division lattice are set to 1/2 chip of code phase), each Doppler frequency shift lattice the chances are the two of KHz Between individual adjacent integers times 2/3 times, code division lattice and a Doppler frequency shift lattice merge becomes one Lattice.Only when the internal carrier wave replicated of receiver and code signal relatively approximate, the output of correlator X_i(τ,Δf_i) just can tend to maximum, and the judgement of the cross-correlation interference of signal checks the power that correlator exports just Reach what maximum realized, if peak power output value exceedes detection thresholding T_NF, then prove that cross-correlation is done Disturb existence.

According to the analysis to GNSS cross-correlated signal, definition GNSS cross-correlation chi square function is:

$X_i(\mathrm{\τ},\mathrm{\Δ}{f}_i)=\frac{1}{N_I}\underset{r=0}{\overset{N_I-1}{\mathrm{\Σ}}}{\left|{\∫}_{r{\mathrm{ZN}}_C}^{(r+1)Z{N}_C}x\left(t\right)c_i^{*}\left(R\left(\mathrm{\Δ}{f}_i\right)(t-\mathrm{\τ})\right)e^{-j2\mathrm{\π\Δ}{f}_{i}t}\mathrm{dt}\right|}^2---\left(16\right)$

Wherein, N_cFor the number of coherent integration chip, N_IFor the number of coherent integration data block, ()^*Represent conjugate complex Number, Ζ is chip-spaced length, and r is integral number of times, and x (t) is required satellite navigation signals, c_i(t) interference letter The chip value of number i, R (Δ f_i) it is that Doppler frequency shift is Δ f under preferable tracking condition_iTime correlation function value.

Due to GNSS signal X₁, X₂...,Obey the sample of normal population N (0,1), then claim statistics Amount

${\mathrm{\χ}}^2=X_1^2+X_2^2+\·\·\·+X_{N_I}^2---\left(17\right)$

Obedience degree of freedom is 2N_Iχ²Distribution.

X_i(τ,Δf_i) be calculated in the grid of a discrete code delay and Doppler frequency shift, if a satellite Correlation peak less than setting thresholding, then cross-correlation interference detection is zero, proceed to next Doppler frequency shift and The grid of code delay composition, if currently detection peak value is more than setting thresholding, then means to exist interference mutually and defends Star signal.

Set up a matrix X_i(n, m), it represents square cross-correlation X_i(τ,Δf_i) value.Code delay τ and Duo Pu Strangle frequency displacement Δ f_iBecome the bivector in matrix.

X_i(n, m)=X_i(τ_n,Δf_im), n=1 ..., N, m=1 ..., M.

Wherein, n represents code phase dimension in a matrix, and m represents Doppler frequency shift dimension in a matrix.

In order to not lose generality, if Doppler frequency shift unit is m₀, at m₀In unit, find and make X_i(n,m₀) Maximum correlation, if n value now is n₀, i.e. n₀=arg max X_i(n,m₀).Finding peak maximum After, with the sampling intervalFind second largest peak value, if n value now isI.e.Then for doppler cells m₀, cross-correlation peak value detector is defined as:

For double alternative hypothesis check problems, observation sample is probably at H₁It is true or at H₀For genuine situation Under draw, it is clear that only both may, according to Bayes's average risk minimum criteria: often obtain a sight Test sample basis, the likelihood ratio of calculating observation sample, i.e.Then with T_NFDetection thresholding enters Row compares, if likelihood ratio is more than or equal to T_NF, it is judged to H₀Otherwise, it is judged to H₁, H₁Represent and disturb letter mutually Number exist situation, H₀Represent the interference non-existent situation of signal, bivector X mutually_i(τ_n,f_dm) it is expressed as 2N_I The χ of degree of freedom²Distribution,For the standard deviation estimate of correlation peak, it is expressed as:

The GNSS cross-correlation peak value detector that the form estimated makes is totally independent of interference noise, has General.In the unit i that code delay and Doppler frequency shift form, if detector income value is less than detection Thresholding, then it is assumed that interfering signal power is little mutually, signal capture is not resulted in impact, directly proceeds to the next one Code delay and the unit i+1 of Doppler frequency shift composition, if detector income value is more than detection thresholding, then recognize Too big with the interfering signal power of Doppler frequency shift unit for code delay, useful signal is caused strong jamming, has Disturb generation mutually, proceed to step 3, disturb elimination mutually, after interference is eliminated, enter next code and prolong Late and Doppler frequency shift component units i+1.

Step 3: determine based on adaptive chip segment signal subspace projection mathematical model.

Assume that GNSS signal y received comprises strong signalAnd weak signalAdditive noise

$\stackrel{\→}{y}=\stackrel{\→}{s}+\stackrel{\→}{w}+\stackrel{\→}{n}---\left(20\right)$

Wherein,For comprising the strong signal of n × 1 vector,For Comprise the weak signal of n × 1 vector,For the matrix of n × m, and m dimension represents m unit amplitude and believes by force Number vector.For the matrix of n × k, and k dimension represents k unit amplitude weak signal vector,For The amplitude vector of the strong signal of m × 1,For the amplitude vector of the weak signal of k × 1,For receiver Noise vector.

The ultimate principle of subspace projection is to utilize the parameter of strong Signal estimation to build strong signal subspace, then obtains To the orthogonal intersection space of its strong signal subspace, then signal will be received on its orthogonal intersection space by projection matrix Project, utilize the orthogonality of the orthogonal intersection space of strong anticipated signal subspace can eliminate strong signal, reach The purpose that cross-correlation interference is suppressed.

According to step 2, when GNSS cross-correlation peak value detector detects interference signal mutually, GNSS connects Receive prow first interference signal mutual to GNSS to capture and follow the tracks of, according to the interference signal parameters obtained, i.e. Doppler frequency shift, code phase and carrier frequency-phase carry out subspace to interference signal and construct and (specifically can be found in 《Fundamentals of Global Positioning System Receivers A Software Approach " book), input signal is projected in the interference space constructed, recycles input signal Deduct the interference signal of estimation, thus obtain the GNSS desired signal needed.

It is P that input signal y projects to disturb signal and required weak signal to construct<HS>subspace mutually_HSY, P_HSy Being further decomposed as two parts, a part projects to<S>subspace disturbing signal to construct mutually, a part Project to the orthocomplement, orthogonal complement subspace for mutual interference spaceAs shown in Figure 2.Disturb Signal estimation mutually Parameter build matrix beThe projection matrix building interference signal mutually isWherein T Represent matrix turns order.

Input signal y projects in the subspace projection that strong signal is constructed, i.e.

$\begin{array}{c}{P}_{s}y=\stackrel{\&LeftRightArrow;}{S}{\left({\stackrel{\&LeftRightArrow;}{S}}^T\stackrel{\&LeftRightArrow;}{S}\right)}^{-1}{\stackrel{\&LeftRightArrow;}{S}}^T\stackrel{\&RightArrow;}{y}\\ =\stackrel{\&LeftRightArrow;}{S}{\left({\stackrel{\&LeftRightArrow;}{S}}^T\stackrel{\&LeftRightArrow;}{S}\right)}^{-1}{\stackrel{\&LeftRightArrow;}{S}}^T(\stackrel{\&LeftRightArrow;}{S}\stackrel{\&LeftRightArrow;}{a_s}+\stackrel{\&LeftRightArrow;}{H}\stackrel{\&RightArrow;}{a_w}+\stackrel{\&RightArrow;}{n})\\ =\stackrel{\&LeftRightArrow;}{S}\stackrel{\&RightArrow;}{a_s}+\stackrel{\&LeftRightArrow;}{S}{\left({\stackrel{\&LeftRightArrow;}{S}}^T\stackrel{\&LeftRightArrow;}{S}\right)}^{-1}{\stackrel{\&LeftRightArrow;}{S}}^T\stackrel{\&LeftRightArrow;}{H}\stackrel{\&LeftRightArrow;}{a_w}+\stackrel{\&LeftRightArrow;}{S}{\left({\stackrel{\&LeftRightArrow;}{S}}^T\stackrel{\&LeftRightArrow;}{S}\right)}^{-1}{\stackrel{\&LeftRightArrow;}{S}}^T\stackrel{\&RightArrow;}{n}\end{array}---\left(21\right)$

Project the strong Signal estimation obtained according to formula (21), required weak signal can be obtained.

$\begin{array}{c}\mathrm{\&Delta;}\stackrel{\&RightArrow;}{y}=\stackrel{\&RightArrow;}{y}-P_s\stackrel{\&RightArrow;}{y}\\ =(I-P_s)\stackrel{\&RightArrow;}{y}={P_s}^{\&perp;}\stackrel{\&RightArrow;}{y}\\ =(\stackrel{\&LeftRightArrow;}{H}-\stackrel{\&LeftRightArrow;}{S}{\left({\stackrel{\&LeftRightArrow;}{S}}^T\stackrel{\&LeftRightArrow;}{S}\right)}^{-1}\left({\stackrel{\&LeftRightArrow;}{S}}^T\stackrel{\&LeftRightArrow;}{H}\right))\stackrel{\&RightArrow;}{a_w}+\mathrm{\&Delta;}\stackrel{\&RightArrow;}{n}\end{array}---\left(22\right)$

The error produced at subspace projection algorithm for receiver noise.

From formula (22), disturb the cross-correlation of signal mutuallyThe property of size direct shadow desired signal Matter, whenTime bigger, each parameter of desired signal can be by serious " distortion ", so that capturing desired signal Time, vacation capture occurs, therefore, uses adaptive subspace projection algorithm to realize the suppression of interference.

Although have good performance based on subspace projection algorithm, but when projection matrix is built Need to carry out cubic matrix multiplying, a matrix inversion.When strong signal number is N, once project fortune When calculating a length of K, altogether need K²(N+1)+2KN²Secondary multiplication K²N+2KN²-NK-N²Sub-addition, and And also need to invert the matrix N × N of, and projection matrix is that a K ties up square formation, needs K²Size Memory space it is stored.Owing to number N of general strong signal is both less than 4, Matrix Calculating now Inverse have formula to use, so matrix inversion is not the subspace projection algorithm difficult point when hardware realizes.Institute With cross-correlation interference Restrainable algorithms based on subspace projection when hardware realizes, if the memory source of system is tight Scarce then need strictly to consider the size of a length of K of each project, if memory source is relatively big, increasing that can be suitable Adding the size of each project length dimension, when making each project, strong signal data is all the most as far as possible In same data bit.

According to the feature of subspace projection algorithm input signal, for convenience of the calculating of subspace projection algorithm, will The data signal of input carries out cutting, is divided into the fragment signal of equal number, uses step 4 conduct the most again Judgment criterion, it is achieved the self-adaptative adjustment of subspace, when output Signal to Interference plus Noise Ratio is more than and presets threshold value, Show that deblocking now is too big, need input data cell time this is finely divided again, it may be assumed that

Wherein, X is whole input signal length, and g is the weighter factor of self adaptation subspace projection, and L is The fragment signal being divided into.

Here, take from the weighter factor adapting to subspace projectionInteger part, K is 1, 2 ... .20, theoretically, if navigation data has saltus step, spread spectrum can be made, output the most no Being continuous wave again, spread spectrum can make capture result be deteriorated.The reason taking 20 is the navigation data at 20ms In, according to GNSS signal pseudo noise code characteristic, at most can only there is a data jump, if front 10ms Data have data jump, then the data of rear 10ms the most do not exist saltus step.Data jump, peaks spectrum Being about 400MHz, such peak value generally can be detected.

Step 4: set up the mutual interference mitigation gains of output chip segment signal, as error judgment criterion again Carry out judging whether to meet criterion, if satisfaction proof obtains desired signal, if being unsatisfactory for returning to Face is rejected again.

The performance of AF panel in step 3 is analyzed, if GNSS desired signal, mutually interference signal and Noise power is respectively ξ_t、ξ_jAnd σ², noise jamming is white Gaussian noise, therefore, ignores noise power, The Signal to Interference plus Noise Ratio receiving signal is:

$R_{\mathrm{in}}=\frac{E\left[{\stackrel{\&RightArrow;}{s}}^H\stackrel{\&RightArrow;}{s}\right]}{E\left[{\stackrel{\&RightArrow;}{w}}^H\stackrel{\&RightArrow;}{w}\right]+E\left[{\stackrel{\&RightArrow;}{n}}^H\stackrel{\&RightArrow;}{n}\right]}=\frac{{\mathrm{\&xi;}}_t}{{\mathrm{\&xi;}}_j+{\mathrm{\&sigma;}}^2}\&ap;\frac{{\mathrm{\&xi;}}_t}{{\mathrm{\&xi;}}_j}---\left(24\right)$

In formula,For strong signal,For weak signal,For additive noise, E [] represents desired signal.

Signal after utilizing self adaptation subspace projection to carry out mutual AF panel is

$\begin{array}{c}\mathrm{\&Delta;}\stackrel{\&RightArrow;}{y}=\stackrel{\&RightArrow;}{y}-P_s\stackrel{\&RightArrow;}{y}\\ =(I-P_s)\stackrel{\&RightArrow;}{y}={P_s}^{\&perp;}\stackrel{\&RightArrow;}{y}\end{array}---\left(25\right)$

Can obtain according to orthogonal intersection cast shadow matrix character, after the mutual AF panel of GNSS, Signal to Interference plus Noise Ratio is

$\begin{array}{c}{R}_{\mathrm{in}}=\frac{E[{\left({P_s}^{\&perp;}\stackrel{\&RightArrow;}{y}\right)}^H\left({p_s}^{\&perp;}\stackrel{\&RightArrow;}{y}\right)}{E\left[{\stackrel{\&RightArrow;}{n}}^H\stackrel{\&RightArrow;}{n}\right]}\\ =\frac{E[{\left({P_s}^{\&perp;}\stackrel{\&RightArrow;}{y}\right)}^H\left({P_s}^{\&perp;}\stackrel{\&RightArrow;}{y}\right)}{{\mathrm{\&sigma;}}^2}\end{array}---\left(26\right)$

Take $\mathrm{\&alpha;}=E\left[{\left({P_s}^{\&perp;}\stackrel{\&RightArrow;}{y}\right)}^H\left({P_s}^{\&perp;}\stackrel{\&RightArrow;}{y}\right)\right]$

Then, 0≤α≤ξ_t, then the mutual interference mitigation gains of GNSS is

G=10lg (R_out/R_in)=10lg [(α ξ_j)/σ²ξ_t] (27)

As α=ξ_tTime, represent that GNSS desired signal is just being fully located at interference signals subspace mutual with GNSS In the subspace handed over, the desired signal power after projection does not lose, and now interference mitigation gains is to the maximum 10lg(^ξ _t/σ²)；When α=0, representing that desired signal is fully located in the mutual interference space of GNSS, son is empty Between project desired signal is produced with to the inhibition as mutually disturbing signal, now, interference mitigation gains It is zero.

Arranging GNSS mutual interference mitigation gains threshold value is G1, due to self adaptation subspace projection weighting because of Son is 20 to the maximum, and 20ms interference remainder error power is generally less than-30dBW, so GNSS disturbs mutually Suppression gain threshold is that G1 takes-30dBW, counts from the weighter factor maximum 20 of self adaptation subspace projection Calculate, draw mutual interference mitigation gains G, when G is less than G1, then each sheet that step (3) is obtained The signal of section combines, and obtains desired signal；If G is more than G1, then return step 3, to certainly The K adapted in the weighter factor of subspace projection carries out add-one operation, to fragment being finely divided again, And return to this step and carry out the calculating of gain, it is thus achieved that new G, then compare with threshold value G1, constantly weigh Multiple, until the G of all of fragment is both less than G1, then the signal of these fragments both less than is combined Obtain desired signal.

The content not being described in detail in description of the invention belongs to the known technology of those skilled in the art.

Claims (3)

1. the mutual disturbance restraining method of self adaptation based on GNSS correlation peak detector, it is characterised in that bag Include following steps:

(1) utilize navigation neceiver to receive aeronautical satellite A and the navigation signal of aeronautical satellite B simultaneously, its The signal of middle aeronautical satellite A is desired signal, and the signal of aeronautical satellite B is interference signal；

(2) between Doppler frequency shift is KHz two adjacent integers times of desired signal and interference signal Code delay and Doppler frequency on set up lattice respectively, and Doppler frequency lattice and a code delay are divided Lattice merge one grid of formation, set up based on χ at each grid²Distribution differentiates hypothesis testing GNSS disturbs correlation peak detector mutually, utilizes GNSS mutual Interference Peaks value detector to carry out cross-correlation interference Detection, without detecting cross-correlation interference, then directly receive desired signal and terminate；If detection To cross-correlation interference, then proceed to step (3)；

Described GNSS disturbs the form of correlation peak detector to be mutually:

$\begin{array}{cc}\frac{X_i(n_0,m_0)-X_i(n_0^o,m_0)}{{\widehat{\mathrm{\&sigma;}}}_{m_0}}\&GreaterEqual;T_{NF}& \begin{array}{c}{H}_0:NoNear-Far\\ H_1:Near-Far\end{array}\end{array}$

Wherein, X_i() is GNSS cross-correlation chi square function,

$X_i(n,m_0)=\frac{1}{N_I}\underset{r=0}{\overset{N_I-1}{\&Sigma;}}{\left|{\&Integral;}_{{\mathrm{rZN}}_C}^{(r+1){\mathrm{ZN}}_C}x\left(t\right)c_i^{*}\left(R\right(m_0\left)\right(t-n\left)\right)e^{-j2{\mathrm{\&pi;m}}_{0}t}dt\right|}^2$

N_CFor the number of coherent integration chip, N_IFor the number of coherent integration data block, ()^*Represent conjugate complex number, Ζ For chip-spaced length, r is integral number of times, and n is code delay, and x (t) is desired signal, c_iT () is interference letter The chip value of number i, R (m₀) it is that Doppler frequency shift is m under preferable tracking condition₀Time correlation function value, m₀ For doppler cells, at m₀In unit, make X_i(n,m₀) n value that maximum correlation is corresponding is n₀, finding After maximum related value, find second largest correlation with sampling interval Δ, make X_i(n,m₀) second largest being correlated with The n value of value correspondence isFor the standard deviation estimate of correlation, T_NFFor detecting thresholding,

${\widehat{\mathrm{\&sigma;}}}_{m_0}=\frac{{\mathrm{\&Sigma;}}_n{X}_i(n,m_0)}{(N-2\mathrm{\&Delta;})\sqrt{N_I}}$

N is the maximum of code phase dimension；

(3) utilize the weighter factor g of self adaptation subspace projection, the data signal received carried out cutting, It is divided into L the fragment signal that length is identical, for each section of fragment signal, does mutually according to what acquisition and tracking arrived Disturb the data bit information of signal, code phase, carrier frequency value and initial phase value, build and disturb letter mutually Number projection matrix, and utilize subspace projection algorithm, from data signal, remove interference signal, obtain each The residual signal of fragment signal enters step (4)；

WhereinX is the total length of the data signal received, and L is the fragment signal number being divided into Amount,Integer part, the value of K be numeral 1～20 in any integer；

(4) calculate the gain of the residual signal of each fragment signal respectively, and compare with gain threshold, right Return in step (3) after the gain fragment signal section of segmentation again beyond gain threshold and utilize subspace Projection algorithm removes interference signal, until the gain of the residual signal of new segmentation section is less than gain threshold；So After by original for the first time without departing from fragment signal and the section of segmentation and again segmenting again of gain threshold After Duan, each segmentation segment signal of the gain respectively less than gain threshold of each segmentation section is combined, it is thus achieved that desired signal.

Self adaptation based on GNSS correlation peak detector the most according to claim 1 is disturbed mutually and is pressed down Method processed, it is characterised in that: the size of described Doppler frequency lattice be Doppler frequency shift be the two of KHz Between individual adjacent integers times 2/3 times.

Self adaptation based on GNSS correlation peak detector the most according to claim 1 and 2 is done mutually Disturb suppressing method, it is characterised in that: 1/2 chip that size is code phase of described code delay lattice.