CN105841714A - High speed X ray pulsar pulse profile delay measurement method - Google Patents

Abstract

The invention discloses a high speed X ray pulsar pulse profile delay measurement method. The method comprises calculating a third-order self-cumulant Csss(t) of a standard pulse profile, a third-order cross-cumulant Csps(t) of a standard pulse profile and an observation pulse profile, calculating a cross-correlation function rxy(tau) of the third-order self-cumulant Csss(t) and the third-order cross-cumulant Csps(t) and estimating pulse time delay through the cross-correlation function. The third-order self-cumulant Csss(t) of the standard pulse profile and the third-order cross-cumulant Csps(t) of the standard pulse profile and the observation pulse profile are calculated so that gaussian noise in a signal is inhibited. The cross-correlation function rxy(tau) replaces a bispectrum so that operation time is reduced and an operation rate is greatly improved. The weighted cross-correlation function spectrum waveform is sharp, a spectrum component of the signal is improved and a precision is obviously improved.

Description

A kind of high speed X-ray pulsar pulse profile time-delay measuring method

Technical field

The present invention relates to Spacecraft Autonomous Navigation Technology field, particularly during a kind of high speed X-ray pulsar pulse profile Prolong measuring method.

Background technology

After first pulsar is found, XNAV (X-ray navigation with independently position) becomes a kind of novel autonomous Airmanship.The X-ray pulsar that widely distributed the cycle is sufficiently stable in immense universe, and spaceborne X penetrates Line photon detector is able to receive that stable pulsar signal.So make XNAV more safe and stable, and have important Construction value and strategic importance, enjoy space flight mechanism of various countries to pay close attention to.U.S. NASA just took the lead in opening a key name early than 2004 " X-ray navigation and autonomous location (XNAV) " project.Europe ESA also opens entitled " ESA deep space probe pulse simultaneously Star navigates " project, carry out numerous studies；In addition, the state such as Russian, Japanese and Chinese has also carried out relevant grinding Study carefully.

TOA determines that (time of advent) basic variable of spacecraft and pulsar position.Whether TOA is estimated accurately straight Connect the performance of impact location.In turn, the precision of TOA is mainly accumulated with calibration pulse star by observation pulsar pile-up pulse profile The precision impact of the time delay between profile.Additionally, time delay can accurately with calibration pulse profile with observe Pulse profile calculates to be estimated.Thus it is guaranteed that the precision of the time delay of X-ray pulsar pulse profile is extremely necessary.

About how improving the precision of pulse TOA or time delay, as Publication about Document proposes some correlation techniques.At document “Emadzadeh A A,Speyer J L.On modeling and pulse phase estimation of X-ray Pulsars.IEEE Transactions on Signal Processing, 2010,58 (9): 4,484 4495 " and document " EA A,Speer J L.X-Ray Pulsar-Based Relative Navigation using Epoch Folding.IEEE Transactions on Aerospace and Electronic Systems, 2011:0018-9251 " middle proposition will arrival Photon carries out folding and obtains folding profile epoch, and the former utilizes Nonlinear least squares fitting, the latter to utilize cross-correlation technique real Now fold the contrast of profile and nominal contour, obtain the value of pulse TOA；Document " Emadzadeh A A, Speyer J L.Asymptotically efficient estimation of pulse time delay for X-ray pulsar Based relative navigation.In:AIAAGN&C Conference, Chicago, IL, 2009:1 12 " middle utilization Maximum likelihood principle, proposes to directly utilize photonic data and tries to achieve the value of pulse TOA by maximizing likelihood function.Rinauro S Et al. at document " Rinauro S, Colonnese S, ScaranoG.Fast near-maximum likelihood phase Estimation of X-ray pulsars.Signal Processing, 2013,93 (1): 326-331 " by pulse TOA in Determine that problem reformulation is a cyclic shift Parameter Estimation Problem, then utilize discrete Fourier transform to complete likelihood function Bigization, obtains the value of pulse TOA；Document " Su Zhe, Xu Luping, Wang Ting etc. a kind of new pulsar pile-up pulse profile time prolongs Late Measurement Algorithm. aerospace journal, 2011,32 (6): 1256-1261 " in propose a kind of rough estimate and combine with accurately measurement XNAV pulse TOA determine method, improve precision by parabola interpolation；Document " thank Zhenhua, Xu Luping, Ni Guangren. Pulsar accumulated pulse profile time delay measurement based on bispectrum. Acta Physica Sinica, 2008,57 (10): 6683-6688 " in profit Suppress Gaussian noise by bispectrum, under low signal-to-noise ratio in the case of also can obtain higher estimated accuracy.Document above design Pulsar TOA method of estimation or impulse time delay method of estimation all do not consider algorithm operation time problems.In spacecraft practical flight In, flight speed is very fast, the important impact estimated on TOA of the length of algorithm operation time.

Therefore, improve algorithm operational efficiency, provided try hard to keep for preferably carrying out precision pulse star Camera calibration Barrier.

Summary of the invention

The technical problem to be solved is, not enough for prior art, it is provided that a kind of high speed X-ray pulsar arteries and veins Rush profile time-delay measuring method.

For solving above-mentioned technical problem, the technical solution adopted in the present invention is:

With three rank from cumulant and the cross-correlation function r of three rank cross cumulant_xy(τ) replace bispectrum to improve computing speed Degree.In the above-mentioned methods, to cross-correlation function r_xy(τ) being weighted, after weighting, the waveform of cross-correlation function becomes more sharp-pointed, letter Spectrum component in number is enhanced, and precision significantly improves.

Compared with prior art, first the had the beneficial effect that present invention of the present invention calculates calibration pulse profile Three rank are from cumulant C_sss(t), calibration pulse profile and three rank cross cumulant C of observation pulse profile_spsT (), thus can press down Gaussian noise in signal processed；The present invention utilizes cross-correlation function r_xy(τ) replace bispectrum, reduce operation time, be greatly improved Arithmetic speed；The present invention is to cross-correlation function r_xy(τ) being weighted, after weighting, the waveform of cross-correlation function becomes more sharp-pointed, letter Spectrum component in number is enhanced, and precision significantly improves.

Accompanying drawing explanation

Fig. 1 is the inventive method flow chart.

Fig. 2 is the comparison diagram of Plays pulse accumulation profile of the present invention and observation pulse accumulation profile；

Fig. 3 is the comparison chart of the present invention and existing method 1 required operation time in an experiment；

Fig. 4 is the present invention and the existing method 1 error condition in the case of different signal to noise ratios, different pulsars；

Detailed description of the invention

The high speed X-ray pulsar pulse profile time-delay measuring method that the present invention provides, can meet X-ray pulsar and lead Boat precision index requirement, thus realize independent navigation.

The inventive method flow diagram as shown in Figure 1, the high speed X-ray pulsar pulse profile time delay of the present invention is surveyed Metering method comprises the following steps:

(1), for selected pulsar, can obtain calibration pulse profile s (n) (wherein n be the burst length/

Impulse phase), and extract cumulative observations pulse profile p (n) by heasoft software, as shown in Figure 2.By following Formula can obtain three rank from cumulant C_sss(t) and three rank cross cumulant C_sps(t):

C_sss(t)=cum [s (n), s (n+t), s (n)]=E [s (n) s (n+t) s (n)] (1)

C_sps(t)=cum [s (n), p (n+t), s (n)]=E [s (n) p (n+t) s (n)]=c_sss(t-τ) (2)

Wherein variable t is observation time, and τ is impulse time delay, and cum is Third-order cumulants operator,

E [s (n) s (n+t) s (n)] is expressed as third moment；

(2), three rank are obtained from cumulant C in step (1) respectively_sss(t) and three rank cross cumulant C_spsT () is carried out from phase Pass computing:

r_xx(τ)=r_xx(t₁-t₂)=r_xx(t₁,t₂)=E [c_sss(t₁)c_sss(t₂)] (3)

r_yy(τ)=r_yy(t₁-t₂)=r_yy(t₁,t₂)=E [c_sps(t₁)c_sps(t₂)] (4)

(3), to three rank in step (1) from cumulant C_sss(t) and three rank cross cumulant C_spsT () carries out computing cross-correlation:

r_xy(τ)=r_xy(t₁-t₂)=r_xy(t₁,t₂)=E [c_sss(t₁)c_sps(t₂)] (5)

(4) peak value of cross-correlation function, can be utilized below to estimate impulse time delay.In step (3) due to background noise Existence, the peak value of cross-correlation function may be affected, and then cause gross error.Therefore, we use SCOT at frequency domain (The Smoothed Coherence Transform) function is to cross-correlation function r_xy(τ) it is weighted, specifically comprises the following steps that

(4a), respectively to r_xy(τ)、r_xx(τ)、r_yy(τ) Fourier transformation is carried out:

$R_{xy}\left(\mathrm{\ω}\right)=\underset{\mathrm{\τ}}{\Σ}{r}_{xy}\left(\mathrm{\τ}\right)e^{-j\mathrm{\ω}\mathrm{\τ}}---\left(6\right)$

$R_{xx}\left(\mathrm{\ω}\right)=\underset{\mathrm{\τ}}{\Σ}{r}_{xx}\left(\mathrm{\τ}\right)e^{-j\mathrm{\ω}\mathrm{\τ}}---\left(7\right)$

$R_{yy}\left(\mathrm{\ω}\right)=\underset{\mathrm{\τ}}{\Σ}{r}_{yy}\left(\mathrm{\τ}\right)e^{-j\mathrm{\ω}\mathrm{\τ}}---\left(8\right)$

Wherein r_xx(τ) it is that standard accumulation profile three rank are from cumulant C_sssThe auto-correlation function of (t), r_yy(τ) it is standard arteries and veins Rush profile and observation pulse profile three rank cross cumulant C_spsThe auto-correlation function of (t), r_xy(τ) it is that three rank are from cumulant C_sss(t) With three rank cross cumulant C_spsT () cross-correlation function, τ is impulse time delay, and ω is angular frequency；

(4b), to the cross-correlation function frequency spectrum R in step (4a)_xy(ω) it is weighted:

$\widehat{R_{xy}\left(\mathrm{\ω}\right)}=R_{xy}\left(\mathrm{\ω}\right)\·H\left(\mathrm{\ω}\right)---\left(9\right)$

Wherein H (ω) is SCOT function (Carter G.C.et al., 1973)；

$H\left(\mathrm{\ω}\right)=\frac{1}{\sqrt{R_{xx}\left(\mathrm{\ω}\right)R_{yy}\left(\mathrm{\ω}\right)}}---\left(10\right)$

R_xx(ω) and R_yy(ω) it is respectively r_xx(τ) and r_yy(τ) frequency spectrum；

(4c), to the R after weighting_xy(ω) inverse Fourier transform is carried out:

$\widehat{r_{xy}\left(\mathrm{\τ}\right)}=\underset{\mathrm{\ω}}{\Σ}\widehat{R_{xy}\left(\mathrm{\ω}\right)}{e}^{j\mathrm{\ω}\mathrm{\τ}}---\left(11\right)$

(4d), finally delay, τ can be obtained by following formula:

$\widehat{\mathrm{\τ}}=\arg \left\{\max \left(\widehat{r_{xy}\left(\mathrm{\τ}\right)}\right)\right\}---\left(12\right)$

Wherein,Represent that searching can makeReach the parameter of peak value

In the present invention, we have selected B0531+21 pulsar as target celestial body.True time delay d=is set 1.0567ms, noise Gaussian distributed.Under this simulated conditions, carry out 100 Monte Carlo Experiments, the present invention can be obtained The average algorithm operation time is 0.697s.

Then utilize document " thank Zhenhua, Xu Luping, Ni Guangren. the pulsar pile-up pulse profile time based on bispectrum prolongs Late measure [J]. Acta Physica Sinica .2008 (10): 6683-6688 " the existing method 1 that proposes, when i.e. utilizing bispectrum to carry out pulse The estimation prolonged, carries out 100 Monte Carlo Experiments equally, and can obtain average operating time is that both 18.772s. specifically contrast Situation is as shown in Figure 3.

By above the simulation result of the inventive method with existing method 1 being contrasted it can be seen that at identical bar Under part, two kinds of methods all can obtain ideal experimental result, but to compare existing method 1 big the operation time of the inventive method Big minimizing, computational efficiency is higher.

Following signal to noise ratio takes-10dB, 0dB, 1dB, 10dB, 100dB respectively, and under the conditions of above different signal to noise ratio values Respectively three pulsars (B0531+21, B1937+21, B1821-24) are carried out 100 Monte Carlo Experiments, obtain the present invention The comparative result of the impulse time delay error of method and existing method 1.Simulation result is it can be seen that the present invention is permissible as shown in Figure 4 Effectively reduce error, improve precision further.

Claims (2)

1. a high speed X-ray pulsar pulse profile time-delay measuring method, it is characterised in that comprise the following steps:

1) three rank of calibration pulse profile are calculated from cumulant C_sssT three rank of (), calibration pulse profile and observation pulse profile are mutual Cumulant C_sps(t)；

2) three rank are calculated from cumulant C_sss(t) and three rank cross cumulant C_spsThe cross-correlation function r of (t)_xy(τ)；

3) to cross-correlation function r_xy(τ) frequency spectrum R_xy(ω) it is weighted obtainingWeight coefficient H (ω) is by r_xx(τ) and r_yy(τ) the spectrum value R after Fourier transformation is carried out_xx(ω) and R_yy(ω) determine；

4) after to weightingCarry out inverse Fourier transform, obtain

5) utilizePeak estimation pulse profile time delay

High speed X-ray pulsar pulse profile time-delay measuring method the most according to claim 1, it is characterised in that described Step 3) the process that implements include:

1) to r_xy(τ)、r_xx(τ)、r_yy(τ) Fourier transformation is carried out:

$R_{xy}\left(\mathrm{\ω}\right)=\underset{\mathrm{\τ}}{\Σ}{r}_{xy}\left(\mathrm{\τ}\right)e^{-j\mathrm{\ω}\mathrm{\τ}}$

$R_{xx}\left(\mathrm{\ω}\right)=\underset{\mathrm{\τ}}{\Σ}{r}_{xx}\left(\mathrm{\τ}\right)e^{-j\mathrm{\ω}\mathrm{\τ}}$

$R_{yy}\left(\mathrm{\ω}\right)=\underset{\mathrm{\τ}}{\Σ}{r}_{yy}\left(\mathrm{\τ}\right)e^{-j\mathrm{\ω}\mathrm{\τ}}$

Wherein r_xx(τ) it is that three rank are from cumulant C_sssThe auto-correlation function of (t), r_yy(τ) it is three rank cross cumulant C_sps(t) from phase Closing function, τ is impulse time delay, and ω is angular frequency；

2) to cross-correlation function r_xy(τ) frequency spectrum R_xy(ω) it is weighted, obtains

$\widehat{R_{xy}\left(\mathrm{\ω}\right)}=R_{xy}\left(\mathrm{\ω}\right)\·H\left(\mathrm{\ω}\right);$

Wherein, H (ω) is SCOT function, i.e. weight coefficient:

$H\left(\mathrm{\ω}\right)=\frac{1}{\sqrt{R_{xx}\left(\mathrm{\ω}\right)R_{yy}\left(\mathrm{\ω}\right)}};$

3) rightCarry out inverse Fourier transform, obtain

$\widehat{r_{xy}\left(\mathrm{\τ}\right)}=\underset{\mathrm{\ω}}{\Σ}\widehat{R_{xy}\left(\mathrm{\ω}\right)}{e}^{j\mathrm{\ω}\mathrm{\τ}}$

4) finally pulse profile time delay can be obtained by following formula

$\widehat{\mathrm{\τ}}=\arg \left\{max\left(\widehat{r_{xy}\left(\mathrm{\τ}\right)}\right)\right\}$

Wherein,Represent that searching can makeReach the parameter of peak value