CN112162152B - Sine wave coherent pulse train signal frequency estimation method based on phase straight line fitting - Google Patents

Sine wave coherent pulse train signal frequency estimation method based on phase straight line fitting Download PDF

Info

Publication number
CN112162152B
CN112162152B CN202010892528.XA CN202010892528A CN112162152B CN 112162152 B CN112162152 B CN 112162152B CN 202010892528 A CN202010892528 A CN 202010892528A CN 112162152 B CN112162152 B CN 112162152B
Authority
CN
China
Prior art keywords
pulse
phase
frequency
value
signal
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202010892528.XA
Other languages
Chinese (zh)
Other versions
CN112162152A (en
Inventor
刘庆云
桂峰
黄启量
王�琦
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Liu Qingyun
Original Assignee
Nanjing Yijieming Information Technology Co ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Nanjing Yijieming Information Technology Co ltd filed Critical Nanjing Yijieming Information Technology Co ltd
Priority to CN202010892528.XA priority Critical patent/CN112162152B/en
Publication of CN112162152A publication Critical patent/CN112162152A/en
Application granted granted Critical
Publication of CN112162152B publication Critical patent/CN112162152B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01RMEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R23/00Arrangements for measuring frequencies; Arrangements for analysing frequency spectra
    • G01R23/02Arrangements for measuring frequency, e.g. pulse repetition rate; Arrangements for measuring period of current or voltage
    • G01R23/12Arrangements for measuring frequency, e.g. pulse repetition rate; Arrangements for measuring period of current or voltage by converting frequency into phase shift
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01RMEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R23/00Arrangements for measuring frequencies; Arrangements for analysing frequency spectra
    • G01R23/16Spectrum analysis; Fourier analysis

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Measuring Phase Differences (AREA)
  • Measuring Frequencies, Analyzing Spectra (AREA)

Abstract

The invention discloses a sine wave coherent pulse train signal frequency estimation method based on phase straight line fitting. The method comprises the following steps: FFT operation is performed on each of a plurality of pulses of an input complex-analysis digital signal, and an average value of initial estimation values of the frequency of each pulse is denoted as f e The method comprises the steps of carrying out a first treatment on the surface of the By f e Performing down-conversion processing on the input complex analysis digital signal, and extracting phase angle information of each pulse; preprocessing the obtained phase angle data and adjusting an initial phase angle; extracting the phase mean value and the center moment of each pulse; performing linear fitting operation on the phase average value of each pulse, and calculating a frequency estimation error delta f according to the result of linear approximation e The method comprises the steps of carrying out a first treatment on the surface of the When Deltaf e When the frequency is not smaller than the set threshold value, the frequency estimated value f e Adjusted to f e =f e +△f e Performing down-conversion again, extracting phase angle information of each pulse again, and repeating the steps; otherwise directly output f e As a final estimate of the frequency. The frequency estimation method has better frequency estimation precision and wider application range.

Description

Sine wave coherent pulse train signal frequency estimation method based on phase straight line fitting
Technical Field
The invention relates to the technical field of digital signal processing, in particular to a sine wave coherent pulse train signal frequency estimation method based on phase straight line fitting.
Background
In many engineering application fields, such as radar/sonar signal processing, communication signal processing, passive positioning, metering, testing, etc., high-precision frequency estimation is required for sine wave coherent pulse train signals. FFT operation is carried out on the signals, and the frequency corresponding to the spectrum peak with the strongest amplitude is taken as the frequencyThe estimated value is a commonly used frequency estimation method (hereinafter, abbreviated as FFT algorithm) of the sine wave coherent pulse train signal. The frequency estimation performance of the FFT algorithm is limited by the digital resolution f of the FFT algorithm in addition to the detection signal-to-noise ratio s /M(f s M is the sampling frequency and the number of FFT operations, respectively). When f s At a certain time, the digital resolution is doubled, and the data storage amount and the calculation amount are correspondingly doubled, so that the FFT algorithm is not suitable for occasions needing high-precision frequency estimation of sine-wave coherent pulse train signals.
To solve this problem, various solutions have been proposed, one of which is: firstly, obtaining a rough frequency estimation value by utilizing an FFT algorithm, and then obtaining a precise frequency estimation value by utilizing a maximum likelihood estimation method within a certain range around the frequency estimation value. When the error of the frequency coarse estimation value is large, the calculation amount of the frequency estimation method is still large. Another effective solution proposed is: firstly, obtaining a rough frequency estimation value by utilizing an FFT algorithm, then, carrying out down-conversion on each coherent pulse by utilizing the estimation value, solving the average value of each pulse signal after down-conversion, then, estimating a frequency error from the average value sequence of each pulse signal by utilizing the FFT algorithm, and correcting the rough frequency estimation value by utilizing the frequency error estimation value (hereinafter, the frequency error estimation value is called as a secondary FFT algorithm). The two-stage FFT algorithm has relatively small calculation complexity, but is only suitable for frequency estimation of coherent pulse trains with fixed heavy frequencies, and the frequency estimation performance is still limited by the detection signal-to-noise ratio and the digital resolution.
Disclosure of Invention
The invention aims to provide a high-precision iteration estimation method for the frequency of a sine wave coherent pulse train signal based on the combination of an FFT algorithm and phase straight line fitting, which can be used for high-precision frequency estimation of the sine wave coherent pulse train signal with arbitrarily changed repetition frequency.
The technical solution for realizing the purpose of the invention is as follows: a sine wave coherent pulse train signal frequency estimation method based on phase straight line fitting comprises the following steps:
step 1, for input complex analysis digitalFFT operation is carried out on a plurality of pulses of the signal, the frequency value corresponding to the strongest spectral line is taken as the primary estimated value of each pulse frequency, and the average value of the primary estimated values of each pulse frequency is taken as f e
Step 2, using the average value f of the primary estimation values of each pulse frequency e Performing down-conversion processing on the input complex analysis digital signal;
step 3, extracting phase angle information of each pulse from the complex analysis digital signal after the down-conversion treatment;
step 4, preprocessing the obtained phase angle data, and adjusting an initial phase angle;
step 5, extracting the phase mean value and the center moment of each pulse;
step 6, performing linear fitting operation on the phase average value of each pulse, and calculating a frequency estimation error delta f according to the result of linear approximation e
Step 7, Δf e Comparing with a set threshold: when Deltaf e When the value is smaller than the set threshold value, directly outputting f e As the final estimated value of the frequency, otherwise, go to step 8;
step 8, frequency estimation value f e Adjusted to f e =f e +Δf e By Δf e And (3) performing down-conversion treatment on the input complex analysis digital signal subjected to down-conversion and initial phase adjustment, and then turning to step (3).
Further, in step 1, the plurality of pulses of the input complex analysis digital signal are respectively subjected to FFT operation, the frequency value corresponding to the strongest spectral line is used as the initial estimation value of each pulse frequency, and the average value of the initial estimation values of each pulse frequency is denoted as f e The method is characterized by comprising the following steps:
recording deviceThe method comprises the steps of analyzing sine wave coherent pulse train signals for noise-containing complex to be processed, wherein T is more than or equal to 0 and less than or equal to T; w (t) is observed white noise with a mean value of zero; f (f) 0 、θ、T、T i 、τ i And P is the signal frequency, initial phase, duration, and center of the ith pulse, respectivelyTime, duration of the ith pulse and number of pulses; />Is of duration tau i The time center is T i Is a rectangular pulse of (2); f (f) 0 The measurement is to be estimated; at the same time, signal s is recorded 0 The sampling frequency and the sampling point number of (t) are f respectively s N, the signal sequence obtained after sampling is s 0 (n),n=1,2,3,…,N-1;
Performing M-point FFT operation on the ith pulse signal, wherein M is greater than the sampling point number of the ith pulse signal, and the frequency value corresponding to the strongest spectral lineAll +.>Then, the average value of the primary estimation values of the pulse frequencies is denoted as f e Then->
Further, step 2 is described as using the average value f of the primary estimation values of each pulse frequency e The down-conversion processing is carried out on the input complex analysis digital signal, and the down-conversion processing is concretely as follows:
by f e For a sequence of sampled signals s 0 (n) down-conversion processing, in which the portions other than the pulse durations are set to zero in the subsequent signal processing, i.eObtain a down-converted sample signal sequence +.>n=1,2,3,…,N-1。
Further, in step 3, the phase angle information of each pulse is extracted from the complex analysis digital signal after the down-conversion processing, which specifically includes the following steps:
computing a data sequencePhase angle information θ of each pulse signal in (a) i (n),n=0,1,2,…,Q i -1, i=1, 2,3, …, P; wherein Q is i The number of sampling points for the ith pulse; the sampling data sequence of the ith pulse is marked as +.>Then->-π≤θ i (n)≤π,n=0,1,2,…,Q i -1。
Further, the step 4 is to preprocess the obtained phase angle data, adjust the initial phase angle, and specifically includes the following steps:
1) Firstly, judging whether the percentage of the phase folding phenomenon occurring in the pulse exceeds a threshold value or not according to the extracted pulse phase data sequences: if the phase folding phenomenon occurs in the pulse phase data of more than 75% of the pulses, the signal initial phase angle is adjusted, the phase data sequence of each pulse is re-extracted, and the specific adjustment method is as follows: order then=0, 1,2, …, N-1; then, extracting the phase data sequence of each pulse; if the pulse number of the phase folding phenomenon of the phase data in the pulse is less than 75%, performing the phase folding operation on the phase data in the pulse with the phase folding phenomenon;
2) The method for eliminating the phase data of which the phase estimation value does not meet the requirement in each pulse comprises the following specific steps: first calculate psi i (n)=θ i (n+1)-θ i (n), n=0, 1,2, …, Q-1, and then solving for the mean value and standard deviation of ψ (n), when the difference between ψ (n) and its mean value is greater than 3 times its standard deviation, either ψ (n-1) or ψ (n+1) is substituted for ψ (n).
Further, the step 5 of extracting the phase mean value and the center moment of each pulse is specifically as follows;
calculating the central moment t of each pulse after phase preprocessing ci And the average value theta of each pulse phase data sequence i Obtaining a new phase data sequence theta i (t ci ),i=1,2,…,P。
Further, in step 6, a line fitting operation is performed on the phase average value of each pulse, and then a frequency estimation error Δf is calculated according to the result of the line fitting e The method is characterized by comprising the following steps:
for phase data sequence theta i (t ci ) Performing linear fitting operation, wherein i=p, p+1, p+2, …, L-1, p and L are index numbers of a start pulse and a final pulse corresponding to a phase data segment with no phase folding at maximum respectively;
assume that the resulting fitted line is Θ (t ci )=2πΔf e t ci0 I=p, p+1, p+2, …, L-1, Δf e I.e. the frequency estimation error.
Further, step 8 is described as Δf e And (3) performing down-conversion treatment on the down-converted and initially-phase-adjusted input complex analysis digital signal again, wherein the formula is as follows:
where n=0, 1,2, …, N-1.
Compared with the prior art, the invention has the remarkable advantages that: (1) The necessary preprocessing operation is carried out on the intra-pulse phase data and inter-pulse phase data, so that the frequency estimation accuracy is better; (2) When phase straight line fitting is performed by utilizing intra-pulse phase data, no constraint condition exists on inter-pulse intervals; (3) The method can be used for high-precision frequency estimation of fixed-repetition-frequency sine-wave coherent pulse train signals, and also can be used for high-precision frequency estimation of random-repetition-frequency sine-wave coherent pulse train signals, and the application range is wider.
Drawings
Fig. 1 is a flowchart of a method for estimating the frequency of a sinusoidal coherent pulse train signal based on phase line fitting according to the present invention.
Fig. 2 is a schematic diagram of a sequence of pulse phase data extracted in step 3 according to an embodiment of the present invention.
Fig. 3 is a schematic diagram of the pulse phase data sequence extracted in step 3 according to the embodiment of the present invention.
Fig. 4 is a schematic diagram of a mean value sequence of the pulse phase data sequences extracted in step 5 according to an embodiment of the present invention.
Fig. 5 is a schematic diagram of the pulse phase data sequence extracted in step 3 after performing a phase line fitting using the phase average sequence shown in fig. 4 according to an embodiment of the present invention.
Fig. 6 is a schematic diagram of a sequence of pulse phase averages extracted from the pulse phase data sequences shown in fig. 4 according to an embodiment of the present invention.
Fig. 7 is a comparison diagram of the estimation performance simulation results of the method and the two-stage FFT algorithm provided by the invention under several input signal-to-noise ratio conditions.
Detailed Description
The invention relates to a sine wave coherent pulse train signal frequency estimation method based on phase straight line fitting, which comprises the following steps:
step 1, FFT operation is respectively carried out on a plurality of pulses of an input complex analysis digital signal, a frequency value corresponding to the strongest spectral line is used as an initial estimation value of each pulse frequency, and an average value of the initial estimation values of each pulse frequency is recorded as f e
Step 2, using the average value f of the primary estimation values of each pulse frequency e Performing down-conversion processing on the input complex analysis digital signal;
step 3, extracting phase angle information of each pulse from the complex analysis digital signal after the down-conversion treatment;
step 4, preprocessing the obtained phase angle data, and adjusting an initial phase angle;
step 5, extracting the phase mean value and the center moment of each pulse;
step 6, performing linear fitting operation on the phase average value of each pulse, and calculating a frequency estimation error delta f according to the result of linear approximation e
Step 7, Δf e Comparing with a set threshold (the threshold is determined by the required estimation accuracy): when Deltaf e When the value is smaller than the set threshold value, directly outputting f e As the final estimated value of the frequency, otherwise, go to step 8;
step 8, frequency estimation value f e Adjusted to f e =f e +Δf e By Δf e And (3) performing down-conversion treatment on the input complex analysis digital signal subjected to down-conversion and initial phase adjustment, and then turning to step (3).
Further, in step 1, the plurality of pulses of the input complex analysis digital signal are respectively subjected to FFT operation, the frequency value corresponding to the strongest spectral line is used as the initial estimation value of each pulse frequency, and the average value of the initial estimation values of each pulse frequency is denoted as f e The method is characterized by comprising the following steps:
recording deviceThe method comprises the steps of analyzing sine wave coherent pulse train signals for noise-containing complex to be processed, wherein T is more than or equal to 0 and less than or equal to T; w (t) is observed white noise with a mean value of zero; f (f) 0 、θ、T、T i 、τ i And P is the signal frequency, initial phase, duration, central time of the ith pulse, duration of the ith pulse and pulse number respectively; i=1, 2,3, …, P;is of duration tau i The time center is T i Is a rectangular pulse of (2); f (f) 0 The measurement is to be estimated; at the same time, signal s is recorded 0 The sampling frequency and the sampling point number of (t) are f respectively s N, the signal sequence obtained after sampling is s 0 (n),n=1,2,3,…,N-1;
Performing M-point FFT operation on the ith pulse signal (M should be greater than the sampling point of the ith pulse signal, but to increase the frequency estimation accuracy as much as possible, and to avoid the inaccuracy of the estimation of each pulse starting and stopping time under the condition of low input signal-to-noise ratio and the introduction of pure receiving noise signal, more conservative pulse starting should be adopted)Time-of-day estimation to ensure that no pure received noise signal portion is introduced), the frequency value corresponding to the strongest spectral lineAll of the pulse frequencies are obtained as initial estimation valuesThen, the average value of the primary estimation values of the pulse frequencies is denoted as f e Then->When the input signal-to-noise ratio is high, for example, greater than 10dB, only the frequency estimation value of a small number of pulses is calculated to reduce the operation amount, and then the average value is obtained.
Further, step 2 is described as using the average value f of the primary estimation values of each pulse frequency e The down-conversion processing is carried out on the input complex analysis digital signal, and the down-conversion processing is concretely as follows:
by f e For a sequence of sampled signals s 0 (n) down-conversion processing, in which the portions other than the pulse durations are set to zero in the subsequent signal processing, i.eObtain a down-converted sample signal sequence +.>n=1,2,3,…,N-1。
Further, in step 3, the phase angle information of each pulse is extracted from the complex analysis digital signal after the down-conversion processing, which specifically includes the following steps:
computing a data sequencePhase angle information θ of each pulse signal in (a) i (n),n=0,1,2,…,Q i -1, i=1, 2,3, …, P; wherein Q is i The number of sampling points for the ith pulse; record the ithThe sampled data sequence of each pulse isThen->-π≤θ i (n)≤π,n=0,1,2,…,Q i -1。
Further, the step 4 is to preprocess the obtained phase angle data, adjust the initial phase angle, and specifically includes the following steps:
1) Firstly, judging whether the percentage of the phase folding phenomenon occurring in the pulse exceeds a threshold value or not according to the extracted pulse phase data sequences: if the phase folding phenomenon occurs in the pulse phase data of more than 75% of the pulses, the signal initial phase angle is adjusted, the phase data sequence of each pulse is re-extracted, and the specific adjustment method is as follows: order then=0, 1,2, …, N-1; then, extracting the phase data sequence of each pulse; if the pulse number of the phase folding phenomenon of the phase data in the pulse is less than 75%, performing the phase folding operation on the phase data in the pulse with the phase folding phenomenon;
2) The method for eliminating the phase data of which the phase estimation value does not meet the requirement in each pulse comprises the following specific steps: first calculate psi i (n)=θ i (n+1)-θ i (n), n=0, 1,2, …, Q-1, and then solving for the mean value and standard deviation of ψ (n), when the difference between ψ (n) and its mean value is greater than 3 times its standard deviation, either ψ (n-1) or ψ (n+1) is substituted for ψ (n).
Further, the step 5 of extracting the phase mean value and the center moment of each pulse is specifically as follows;
calculating the central moment t of each pulse after phase preprocessing ci And the average value theta of each pulse phase data sequence i Obtaining a new phase data sequence theta i (t ci ),i=1,2,…,P。
Further, in step 6, a straight line fitting operation is performed on the phase average value of each pulse, and thenCalculating frequency estimation error Deltaf according to the result of straight line similarity e The method is characterized by comprising the following steps:
for phase data sequence theta i (t ci ) Performing linear fitting operation, wherein i=p, p+1, p+2, …, L-1, p and L are index numbers of a start pulse and a final pulse corresponding to a phase data segment with no phase folding at maximum respectively;
assume that the resulting fitted line is Θ (t ci )=2πΔf e t ci0 I=p, p+1, p+2, …, L-1, Δf e I.e. the frequency estimation error.
The phase data sequence theta obtained in the step 6 can also be used for the following steps before the step 6 is carried out i (t ci ) (i=1, 2, …, P) performs a linearization operation, i.e., an unfolding operation, so that more pulse phase data participates in the line fitting operation.
Further, step 8 is described as Δf e And (3) performing down-conversion treatment on the down-converted and initially-phase-adjusted input complex analysis digital signal again, wherein the formula is as follows:
where n=0, 1,2, …, N-1.
Compared with a two-stage FFT algorithm, the frequency estimation accuracy is better although the calculation complexity is increased, and the limitation of the two-stage FFT algorithm is overcome, so that the method can be used for high-accuracy frequency estimation of a fixed-repetition-frequency sine-wave coherent pulse train signal and high-accuracy frequency estimation of a sine-wave coherent pulse train signal with random change of the repetition frequency, and the application range is wider.
The invention is described in further detail below with reference to the accompanying drawings and specific examples.
Examples
Fig. 1 shows specific processing steps of a single sine wave coherent burst signal frequency estimation method employed in the present invention. The noise-containing complex analysis sine type coherent pulse train signal to be processed in the invention can be expressed as a form shown in (1).
In the formula (1), w (t) is observed white noise with a mean value of zero; f (f) 0 、θ、T、T i 、τ i And P is the signal frequency, initial phase, duration, central time of the ith pulse, duration of the ith pulse and pulse number respectively;is of duration tau i The time center is T i Is a rectangular pulse of (2); f (f) 0 Is the quantity to be estimated. At the same time, signal s is recorded 0 The sampling frequency and the sampling point number of (t) are f respectively s N, the signal sequence obtained after sampling is s 0 (N) (n=1, 2,3, …, N-1). The phase of the signal sequence can be expressed in the form shown in formula (2) (although the invention processes the sampled digital signal, the correspondence between signal phase and time is represented by the continuous time variable t, which does not affect the description of the invention).
In the formula (2), the amino acid sequence of the compound,is the phase noise generated by w (t). When w (t) is random noise with mean value zero, ++>Also random noise with zero mean. Theoretically, ignore +.>Phi (t) varies linearly with time t, both during the respective pulse duration and during the whole observation time. But in the opposite directionWhen phi (t) is calculated numerically, phi (t) can only take values within the (-pi, pi) range, so that phi (t) is in fact periodically and linearly changed in the whole observation time unless the sine wave coherent pulse train signal to be processed is changed in only one period in the whole observation time. When->Cannot be ignored unless +.>Is strong enough to drown out phi (t), otherwise such a law of variation of phi (t) should still be observed.
Performing M-point FFT operation on the ith (i=1, 2,3, …, P) pulse signal (M should be greater than the sampling point of the ith pulse signal, but should increase the value of M as much as possible to improve the frequency estimation accuracy; to avoid the inaccurate estimation of the starting and ending moments of each pulse under the condition of lower input signal-to-noise ratio and the introduction of pure receiving noise signals, more conservative estimation of the starting and ending moments of the pulse should be adopted to ensure that the part of the pure receiving noise signals is not introduced), obtaining the frequency value corresponding to the strongest spectral line, and recording asAll are obtainedThen, obtain the coarse frequency estimation value +.>
By f e For signal sequence s 0 (N) (n=1, 2,3, …, N-1) after the down-conversion processing, the extracted signal phase information can be expressed in the form of expression (3).
In formula (3), Δf=f 0 -f e . Due to frequency estimationThe gauge error Δf is typically much smaller than f 0 So, the extracted signal phase information may no longer be periodically varying, at least the period of variation will increase significantly.
Subject to phase noiseAnd the influence of the initial value theta of the phase on the signal sequence s 0 (N) (n=1, 2,3, …, N-1) there is a possibility that the intra-pulse phase information extracted after the down-conversion processing has out-of-tolerance or phase-folded data. The necessary preprocessing of the extracted phase data is a necessary condition for obtaining a high-precision frequency estimate.
FIG. 2 shows the embodiment in using f e For signal sequence s 0 (N) (n=1, 2,3, …, N-1) after the down-conversion processing (the down-converted signal sequence is recorded) the extracted pulse phase data sequences. As can be seen from fig. 2, the extracted pulse phase data sequences all have a phase folding phenomenon. In this case, the down-converted signal sequence is multiplied by a factor e Then, the phase data sequence of each pulse is re-extracted. FIG. 3 shows the embodiment in using f e For signal sequence s 0 (N) (n=1, 2,3, …, N-1) each pulse phase data sequence extracted after the down-conversion processing. As can be seen from fig. 3, the phase data sequences of the second half pulse all have phase folding phenomenon, while the phase data sequences of the first half pulse have no phase folding phenomenon, but have a problem that the individual phase estimation values are significantly out of tolerance. In this case, for the problem that only the individual phase estimates are significantly out of tolerance within the pulse, the estimate may be eliminated from participating in the subsequent calculation process. For the data with phase folding and individual out-of-tolerance in pulse, the phase folding should be firstly released, and then the individual out-of-tolerance data should be removed.
After preprocessing the extracted intra-pulse phase data, calculating an intra-pulse phase average value of each pulse to obtain an intra-pulse phase average value data sequence
In the formula (4), the amino acid sequence of the compound,is the noise part in the phase mean, +.>FIG. 4 is an intra-pulse phase average sequence extracted after preprocessing the intra-pulse phase data sequence shown in FIG. 3 according to the present invention. As can be seen from fig. 4, this phase-averaged sequence has a phase folding phenomenon.
In straight line fitting the pulse phase averaging sequence shown in fig. 4, two methods can be employed: firstly, performing an unfolding phase folding operation on an intra-pulse phase mean value sequence with a phase folding phenomenon, and then performing straight line fitting; and secondly, taking the intra-pulse phase average value sequence with the maximum phase folding, and performing straight line fitting on the first 11 intra-pulse phase average value sequences in fig. 4.
Performing straight line fitting on the intra-pulse phase average value sequence to obtain an estimated value delta f of the frequency difference e After that, if Δf e Less than or equal to a given threshold value, the frequency estimate f is output e The frequency estimation process ends. When Deltaf e When the frequency is greater than a given threshold value, the frequency estimation value is adjusted according to the formula (5), and delta f is utilized e For down-converted signal sequence(n=0, 1,2, …, N-1) performing down-conversion again, and then sequentially performing extraction, pretreatment, intra-pulse phase average value sequence extraction and frequency difference re-iteration estimation on each pulse phase data sequence until the frequency difference estimated value is smaller than a given threshold. Generally, the frequency estimation process can be ended through two iterative operations of the frequency error estimation value.
Fig. 5 and 6 are graphs showing the frequency offset estimation Δf obtained by performing a phase straight line fitting using the phase average data of the first 11 pulses shown in fig. 4 e (accordingly, f e =f e +Δf e ) The estimated value is used for down-converting the signal sequence(n=0, 1,2, …, N-1) the pulse phase data sequences extracted after the re-down conversion and intra-pulse phase mean value sequences.
FIG. 7 shows comparison diagrams of estimation performance simulation results of the method and the secondary FFT algorithm under several input signal-to-noise ratio conditions (signal sampling frequency, coherent pulse train duration, pulse repetition frequency and pulse width are respectively 100MHz, 1ms, 20kHz and 2 μs; the frequency of a signal to be estimated is uniformly valued within the range of 5 MHz+/-1.05 kHz; the number of simulation times under the condition of single input signal-to-noise ratio is 500), the same frequency coarse estimation value is used in the method and the secondary FFT algorithm, the digital resolution of the FFT algorithm used is 7.3Hz when the frequency error estimation is carried out in the secondary FFT algorithm, and the starting and ending time of each pulse is assumed to be known when the simulation calculation is carried out. As can be seen from fig. 7, in a given input signal-to-noise ratio range, the estimation performance of the frequency estimation method provided by the present invention is better than that of the two-stage FFT algorithm; as the input signal-to-noise ratio decreases, the estimation performance of the two frequency estimation methods gets closer and closer.
The invention carries out necessary preprocessing operation on the intra-pulse and inter-pulse phase data, and has better frequency estimation precision; when phase straight line fitting is performed by utilizing intra-pulse phase data, no constraint condition exists on inter-pulse intervals; the method can be used for high-precision frequency estimation of fixed-repetition-frequency sine-wave coherent pulse train signals, and also can be used for high-precision frequency estimation of random-repetition-frequency sine-wave coherent pulse train signals, and the application range is wider.

Claims (7)

1. A sine wave coherent pulse train signal frequency estimation method based on phase straight line fitting is characterized by comprising the following steps:
step 1, opposite transfusionPerforming FFT operation on multiple pulses of the input complex analysis digital signal, taking the frequency value corresponding to the strongest spectral line as the initial estimation value of each pulse frequency, and recording the average value of the initial estimation values of each pulse frequency as f e
Step 2, using the average value f of the primary estimation values of each pulse frequency e Performing down-conversion processing on the input complex analysis digital signal;
step 3, extracting phase angle information of each pulse from the complex analysis digital signal after the down-conversion treatment;
and 4, preprocessing the obtained phase angle data, and adjusting an initial phase angle, wherein the method comprises the following steps of:
1) Firstly, judging whether the percentage of the phase folding phenomenon occurring in the pulse exceeds a threshold value or not according to the extracted pulse phase data sequences: if the phase folding phenomenon occurs in the pulse phase data of more than 75% of the pulses, the signal initial phase angle is adjusted, the phase data sequence of each pulse is re-extracted, and the specific adjustment method is as follows: down-converted sample signal sequenceThen, extracting the phase data sequence of each pulse; if the pulse number of the phase folding phenomenon of the phase data in the pulse is less than 75%, performing the phase folding operation on the phase data in the pulse with the phase folding phenomenon;
2) The method for eliminating the phase data of which the phase estimation value does not meet the requirement in each pulse comprises the following specific steps: first calculate psi i (n)=θ i (n+1)-θ i (n), n=0, 1,2, …, Q-1, i=1, 2,3, …, P representing the number of pulses, θ i (n) is a sequencePhase angle information of each pulse signal; then solving the mean value and standard deviation of the phi (n), and replacing the phi (n) by the phi (n-1) or the phi (n+1) when the difference between the phi (n) and the mean value is more than 3 times of the standard deviation;
step 5, extracting the phase mean value and the center moment of each pulse;
step (a)6, carrying out linear fitting operation on the phase average value of each pulse, and calculating a frequency estimation error delta f according to the linear fitting result e
Step 7, Δf e Comparing with a set threshold: when Deltaf e When the value is smaller than the set threshold value, directly outputting f e As the final estimated value of the frequency, otherwise, go to step 8;
step 8, frequency estimation value f e Adjusted to f e =f e +Δf e By Δf e And (3) performing down-conversion treatment on the input complex analysis digital signal subjected to down-conversion and initial phase adjustment, and then turning to step (3).
2. The method for estimating frequency of sine wave coherent pulse train signal based on phase straight line fitting as set forth in claim 1, wherein in step 1, the plurality of pulses of the input complex analysis digital signal are subjected to FFT operation, and a frequency value corresponding to a strongest spectral line is used as an initial estimate value of each pulse frequency, and an average value of the initial estimate values of each pulse frequency is denoted as f e The method is characterized by comprising the following steps:
recording deviceThe method comprises the steps of analyzing sine wave coherent pulse train signals for noise-containing complex to be processed, wherein T is more than or equal to 0 and less than or equal to T; w (t) is observed white noise with a mean value of zero; f (f) 0 、θ、T、T i 、τ i And P is the signal frequency, initial phase, duration, central time of the ith pulse, duration of the ith pulse and pulse number respectively; i=1, 2,3, …, P; />Is of duration tau i The time center is T i Is a rectangular pulse of (2); f (f) 0 The measurement is to be estimated; at the same time, signal s is recorded 0 The sampling frequency and the sampling point number of (t) are f respectively s N, the signal sequence obtained after sampling is s 0 (n),n=1,2,3,…,N-1;
For the ith pulse signalM-point FFT operation, M is greater than the sampling point number of the ith pulse signal, and the frequency value corresponding to the strongest spectral lineAll +.>Then, the average value of the primary estimation values of the pulse frequencies is denoted as f e Then->
3. The method for estimating a frequency of a sinusoidal coherent pulse train signal based on phase straight line fitting according to claim 2, wherein the average value f of the initial estimation values of each pulse frequency is used in step 2 e The down-conversion processing is carried out on the input complex analysis digital signal, and the down-conversion processing is concretely as follows:
by f e For a sequence of sampled signals s 0 (n) down-conversion processing, in which the portions other than the pulse durations are set to zero in the subsequent signal processing, i.eObtain a down-converted sample signal sequence +.>
4. The method for estimating frequency of sinusoidal coherent pulse train signal based on phase straight line fitting according to claim 3, wherein in step 3, the phase angle information of each pulse is extracted from the complex analysis digital signal after the down-conversion processing, specifically as follows:
computing a data sequencePhase angle information θ of each pulse signal in (a) i (n),n=0,1,2,…,Q i -1, i=1, 2,3, …, P; wherein Q is i The number of sampling points for the ith pulse; the sampling data sequence of the ith pulse is marked as +.>Then->-π≤θ i (n)≤π,n=0,1,2,…,Q i -1。
5. The method for estimating a frequency of a sinusoidal coherent pulse train signal based on phase straight line fitting according to claim 4, wherein the step 5 is characterized by extracting a phase mean value and a center moment of each pulse, specifically as follows;
calculating the central moment t of each pulse after phase preprocessing ci And the average value theta of each pulse phase data sequence i Obtaining a new phase data sequence theta i (t ci ),i=1,2,…,P。
6. The method for estimating a frequency of a sinusoidal coherent pulse train signal based on phase straight line fitting according to claim 5, wherein in step 6, a straight line fitting operation is performed on the phase average value of each pulse, and then a frequency estimation error Δf is calculated according to the result of the straight line fitting e The method is characterized by comprising the following steps:
for phase data sequence theta i (t ci ) Performing linear fitting operation, wherein i=p, p+1, p+2, …, L-1, p and L are index numbers of a start pulse and a final pulse corresponding to a phase data segment with no phase folding at maximum respectively;
assume that the resulting fitted line is Θ (t ci )=2πΔf e t ci0 I=p, p+1, p+2, …, L-1, Δf e I.e. the frequency estimation error.
7. According to claimThe method for estimating a frequency of a sinusoidal coherent pulse train signal based on phase line fitting as set forth in claim 6, wherein said step 8 uses Δf e And (3) performing down-conversion treatment on the down-converted and initially-phase-adjusted input complex analysis digital signal again, wherein the formula is as follows:
where n=0, 1,2, …, N-1.
CN202010892528.XA 2020-08-31 2020-08-31 Sine wave coherent pulse train signal frequency estimation method based on phase straight line fitting Active CN112162152B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010892528.XA CN112162152B (en) 2020-08-31 2020-08-31 Sine wave coherent pulse train signal frequency estimation method based on phase straight line fitting

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010892528.XA CN112162152B (en) 2020-08-31 2020-08-31 Sine wave coherent pulse train signal frequency estimation method based on phase straight line fitting

Publications (2)

Publication Number Publication Date
CN112162152A CN112162152A (en) 2021-01-01
CN112162152B true CN112162152B (en) 2024-01-26

Family

ID=73859414

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010892528.XA Active CN112162152B (en) 2020-08-31 2020-08-31 Sine wave coherent pulse train signal frequency estimation method based on phase straight line fitting

Country Status (1)

Country Link
CN (1) CN112162152B (en)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113093123B (en) * 2021-04-06 2023-08-08 南京工程学院 Jammer for resisting pulse Doppler radar and interference method thereof
CN113406386B (en) * 2021-06-23 2023-04-25 中国电子科技集团公司第二十九研究所 Signal frequency accurate estimation method based on digital down-conversion
CN115856424A (en) * 2023-03-01 2023-03-28 西安瀚博电子科技有限公司 Signal frequency and amplitude self-adaptive extraction method based on peak-to-adjacent ratio

Citations (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6839381B1 (en) * 2000-01-12 2005-01-04 Freescale Semiconductor, Inc. Method and apparatus for coherent detection in a telecommunications system
CN101984613A (en) * 2010-11-30 2011-03-09 中国工程物理研究院电子工程研究所 Code rate estimation method for low-rate binary phase shift keying (BPSK) burst signal
CN103188067A (en) * 2013-04-03 2013-07-03 西安星河亮点信息技术有限公司 Method for estimating and correcting deviation and error of chip clock frequency of spread spectrum system
CN105302935A (en) * 2015-08-10 2016-02-03 工业和信息化部电信研究院 Digital demodulating and measurement analysis method
CN106872773A (en) * 2017-04-25 2017-06-20 中国电子科技集团公司第二十九研究所 A kind of the multiple-pulse Precision Method of Freuqency Measurement and device of single carrier frequency pulse signal
CN109061296A (en) * 2018-07-17 2018-12-21 南京恒电电子有限公司 A kind of high-precision carrier frequency estimation method of RF pulse signal
CN109270344A (en) * 2018-10-07 2019-01-25 扬州大学 Coherent pulse signal frequency estimating methods under pulse missing
CN109342813A (en) * 2018-12-24 2019-02-15 常州工学院 A kind of sinusoidal signal frequency estimation method based on DFT and dichotomy
CN109738697A (en) * 2019-01-24 2019-05-10 中国电子科技集团公司第二十九研究所 A kind of frequency measurement method based on the correction of finite point discrete spectrum

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
GB0103669D0 (en) * 2001-02-15 2001-03-28 Central Research Lab Ltd A method of estimating the carrier frequency of a phase-modulated signal
US9450598B2 (en) * 2015-01-26 2016-09-20 Guzik Technical Enterprises Two-stage digital down-conversion of RF pulses

Patent Citations (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6839381B1 (en) * 2000-01-12 2005-01-04 Freescale Semiconductor, Inc. Method and apparatus for coherent detection in a telecommunications system
CN101984613A (en) * 2010-11-30 2011-03-09 中国工程物理研究院电子工程研究所 Code rate estimation method for low-rate binary phase shift keying (BPSK) burst signal
CN103188067A (en) * 2013-04-03 2013-07-03 西安星河亮点信息技术有限公司 Method for estimating and correcting deviation and error of chip clock frequency of spread spectrum system
CN105302935A (en) * 2015-08-10 2016-02-03 工业和信息化部电信研究院 Digital demodulating and measurement analysis method
CN106872773A (en) * 2017-04-25 2017-06-20 中国电子科技集团公司第二十九研究所 A kind of the multiple-pulse Precision Method of Freuqency Measurement and device of single carrier frequency pulse signal
CN109061296A (en) * 2018-07-17 2018-12-21 南京恒电电子有限公司 A kind of high-precision carrier frequency estimation method of RF pulse signal
CN109270344A (en) * 2018-10-07 2019-01-25 扬州大学 Coherent pulse signal frequency estimating methods under pulse missing
CN109342813A (en) * 2018-12-24 2019-02-15 常州工学院 A kind of sinusoidal signal frequency estimation method based on DFT and dichotomy
CN109738697A (en) * 2019-01-24 2019-05-10 中国电子科技集团公司第二十九研究所 A kind of frequency measurement method based on the correction of finite point discrete spectrum

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
相参脉冲串多普勒频率变化率估计算法;张刚兵 等;数据采集与处理;第26卷(第2期);第135-139页 *
雷达信号相参脉冲串频率测量方法研究;卢鑫 等;航天电子对抗;第25卷(第6期);第30-32、61页 *

Also Published As

Publication number Publication date
CN112162152A (en) 2021-01-01

Similar Documents

Publication Publication Date Title
CN112162152B (en) Sine wave coherent pulse train signal frequency estimation method based on phase straight line fitting
CN109471095B (en) FMCW radar distance estimation method based on fast iterative interpolation
AU2010340812B2 (en) Pulse radar range profile motion compensation
CN109490862B (en) Carrier frequency estimation method based on phase difference statistical spectrum
CN111562438B (en) Sinusoidal signal frequency estimation method and device based on FFT and phase difference
CN112881796A (en) Multi-frequency real signal frequency estimation algorithm for spectrum leakage correction
CN109682492B (en) Frequency estimation method based on frequency domain Gaussian fitting
CN111865865A (en) Frequency offset and phase offset estimation method suitable for high-sensitivity satellite-borne ADS-B receiver
CN114895248A (en) Sinusoidal frequency modulation signal parameter estimation method, system and medium
CN113406386B (en) Signal frequency accurate estimation method based on digital down-conversion
CN110441749B (en) Frequency stepping radar target motion parameter estimation method
CN112883318A (en) Multi-frequency attenuation signal parameter estimation algorithm of subtraction strategy
CN111007473A (en) High-speed weak target detection method based on distance frequency domain autocorrelation function
CN106777505A (en) The frequency estimating methods and device of the robust of the undersampled signal based on frequency deviation identification
CN110808929A (en) Real-complex conversion type signal-to-noise ratio estimation algorithm of subtraction strategy
CN112964931B (en) Non-ideal multi-damping harmonic signal parameter measurement method based on two-channel undersampling
CN112505413B (en) Time-frequency analysis method and system
CN112162153B (en) Sine wave signal frequency estimation method based on phase straight line fitting
CN112129983B (en) Waveform recovery data processing method based on equivalent sampling at equal time intervals
CN111487593B (en) Incomplete radar signal restoration method
CN112748285A (en) Phase measurement method based on intelligent tracking correlation operation
CN108132460B (en) Pulse compression compensation algorithm based on frequency domain channel equalization
CN103441975A (en) Two-phase coding signal parameter estimation method based on power spectrum
CN107576842B (en) Broadband synchronous sampling method
CN116449304B (en) SAR emission pulse arrival time measurement method based on frequency measurement

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
TR01 Transfer of patent right

Effective date of registration: 20240415

Address after: Building 10, Love Garden, No. 6 Hongtai Road, Pukou District, Nanjing City, Jiangsu Province, 210000

Patentee after: Liu Qingyun

Country or region after: China

Address before: 210039 room 1216, 12 / F, building 33, Chaoyang Xiyuan business building, Banqiao street, Yuhuatai District, Nanjing City, Jiangsu Province

Patentee before: Nanjing yijieming Information Technology Co.,Ltd.

Country or region before: China

TR01 Transfer of patent right