CN101626357B - Carrier synchronization method of MPSK system based on maximum likelihood estimation - Google Patents

Abstract

The invention discloses a carrier synchronization method of an MPSK system based on maximum likelihood estimation, belonging to the technical field of the modulation and the demodulation of digital communication and comprising the following steps: firstly, eliminating the modulation information of a receiving signal to obtain a single carrier only polluted by noises; then, utilizing a fast search algorithm based on dichotomy to obtain maximum likelihood estimation values of frequency deviation and phase deviation; and then, compensating an original receiving signal by the maximum likelihood estimation values to complete carrier synchronization and demodulation, wherein the estimation and compensation method needs to be continuously used for tracking the carrier synchronization in the demodulation process. A simulation result shows that the carrier synchronization method has very good property approximated to the Cramer-Rao lower bound, very little calculation amount and very wide frequency deviation estimation range and is a practical carrier synchronization method.

Description

A kind of MPSK system carrier method for synchronous based on maximal possibility estimation

Technical field

The present invention relates to the carrier synchronization method of M level phase shift keying (MPSK) system, belong to the modulation-demodulation technique field in the digital communication.

Background technology

In satellite communication, deep space communication and various wireless communication system, mpsk signal is modulated because of having permanent envelope, power validity advantages of higher, is one of the most frequently used modulation system.Because the employing of various advanced chnnel codings makes communication system under low signal-to-noise ratio, also have good performance, but also the MPSK system is had higher requirement simultaneously, promptly can reliable work under low signal-to-noise ratio.Wherein, needing the problem of the most critical of solution is the carrier synchronization problem.

Communication system comprises two kinds of transmission modes: burst communication and continuous communiction.The data of burst communication are made up of frame head and data, and frame head is used for the beginning of mark burst frame with synchronously, and its content is known for receiving terminal, and the burst frame data format is as shown in Figure 1.The data of continuous communiction then are to send continuously; Middle duties ground is sent out some and is used to remove the unique word of phase ambiguity; The data content of unique word also is known for receiving terminal; But general unique word is very short, under the low signal-to-noise ratio condition, can't be used for carrier synchronization, and the continuous communiction data format is as shown in Figure 2.

For realizing the carrier synchronization of mpsk signal, at present the most frequently used is PHASE-LOCKED LOOP PLL TECHNIQUE, i.e. section's Stas (Costas) ring.Be modulated to example with BPSK; The schematic diagram of section's Stas ring is as shown in Figure 3, and in this loop, voltage controlled oscillator (VCO) provides two-way mutually orthogonal carrier wave; In homophase and two phase discriminators of quadrature, carry out phase demodulation respectively with the bpsk signal of input, behind low pass filter, obtain signal v ₅, v ₆, deliver to a multiplier again and multiply each other, remove v ₅, v ₆In digital signal, obtain reflecting the error controling signal v of the difference of VCO and incoming carrier phase place ₇Suppose that loop is locked, if consideration of noise not, then the input signal r (t) of loop does

r(t)＝x(t)cosω _ct (1)

The local reference signal v of homophase and two phase discriminators of quadrature then ₁, v ₂Be respectively

$\left\{\begin{array}{c}{v}_1=\cos ({\mathrm{\ω}}_{c}t+\mathrm{\θ})\\ v_2=\sin ({\mathrm{\ω}}_{c}t+\mathrm{\θ})\end{array}\right.---\left(2\right)$

Input signal r (t) and v ₁, v ₂After multiplying each other respectively, correspondence obtains the signal v behind the phase demodulation ₃And v ₄:

$\left\{\begin{array}{c}{v}_3=x\left(t\right)\cos {\mathrm{\ω}}_{c}t\cos ({\mathrm{\ω}}_{c}t+\mathrm{\θ})=\frac{1}{2}x\left(t\right)[\cos \mathrm{\θ}+\cos (2{\mathrm{\ω}}_{c}t+\mathrm{\θ})]\\ v_4=x\left(t\right)\cos {\mathrm{\ω}}_{c}t\sin ({\mathrm{\ω}}_{c}t+\mathrm{\θ})=\frac{1}{2}x\left(t\right)[\sin \mathrm{\θ}+\sin (2{\mathrm{\ω}}_{c}t+\mathrm{\θ})]\end{array}\right.---\left(3\right)$

v ₃, v ₄Respectively through getting behind the low pass filter

$\left\{\begin{array}{c}{v}_5=\frac{1}{2}x\left(t\right)\cos \mathrm{\θ}\\ v_6=\frac{1}{2}x\left(t\right)\sin \mathrm{\θ}\end{array}\right.---\left(4\right)$

v ₅, v ₆After the process multiplier multiplies each other,

$v_7=v_5{v}_6=\frac{1}{4}{x}^2\left(t\right)\sin \mathrm{\θ}\cos \mathrm{\θ}=\frac{1}{8}{x}^2\left(t\right)\sin 2\mathrm{\θ}\≈\frac{1}{4}{x}^2\left(t\right)\mathrm{\θ}---\left(5\right)$

This voltage is controlled VCO later on through loop filter, makes it and ω _cWith frequently, phase place only differs from a very little θ.This moment v ₁=cos (ω _cT+ θ) be the sync carrier that will extract, and ${\mathrm{<figure-callout\; id="5"\; label="v"\; filenames="HA20178328200910087839201E00032.png"\; state="\{\{state\}\}">v</figure-callout>}}_5=\frac{1}{2}x\left(t\right)\mathrm{Cos}\mathrm{\&theta;}\&ap;\frac{1}{2}x\left(t\right)$ Output for demodulator.

Can find out from said process, work as v ₁=cos (ω _cDuring t+ θ+π), though the formula that obtains (5) has identical result, the demodulation result that obtains has been inverted fully, and promptly there is phase ambiguity in section's Stas ring.For BPSK 2 ambiguityes are arranged, M ambiguity then arranged for MPSK.

In actual applications, fuzziness can be eliminated by unique word, because unique word is known for receiving terminal, by unique word is compared with demodulation result, just can remove phase ambiguity.But, under the low signal-to-noise ratio condition, section's Stas ring to go into the lock time long, be not suitable for burst communication.And when signal to noise ratio was lower than 6dB, the cycle-skipping phenomenon of section's Stas ring in tracing process (phase ambiguity appears again in the carrier wave of promptly removing after the phase ambiguity) was just apparent in view, no longer is suitable in the coherent demodulator.

In order to make full use of the information of data, improve the performance of carrier synchronization, at present, and more existing method for synchronous based on frequency offset estimating and compensation, it realizes that basically principle is as shown in Figure 4.The intermediate-freuqncy signal that receives at first changes to base band through the local oscillator quadrature frequency conversion, passes through the AD sample varianceization again.Supposing has desirable bit timing, then the baseband signal r after down-conversion _kCan be expressed as

r _k＝exp(jθ _k)exp(j2πf _dkT _s+jφ ₀)+z _k (6)

In the formula (6), θ _kBe the phase place of constellation point symbol of transmission, i.e. modulation intelligence, for MPSK, its value is θ _k∈ 0,2 π/M ..., 2 (M-1) π/M}, M is an order of modulation; f _dIt is residual frequency departure; T _sThe is-symbol cycle; φ ₀It is initial skew; z _kBe noise, be without loss of generality that suppose that it is additivity white complex gaussian noise (AWGN), its average is 0, variance is σ ²J is an imaginary unit; K representes the label of k sampled point.The power normalization of signal is 1.

Obviously, as long as with residual frequency departure f _dWith initial skew φ ₀Estimation is come out, and its estimated value is expressed as respectively

With

Then can recover original signal, obtain demodulation result d _k:

$d_k=r_k\exp (-\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;}{\hat{f}}_{d}k{T}_s-j{\widehat{\mathrm{\&phi;}}}_0)---\left(7\right)$

But these methods of estimation, perhaps operand is too big, and perhaps estimated accuracy is not high enough, is not suitable for the situation of low signal-to-noise ratio, and the frequency deviation region that perhaps adapts to is limited, all is not suitable for Project Realization.

Summary of the invention

The objective of the invention is deficiency, propose a kind of MKSP system carrier method for synchronous based on maximal possibility estimation to existing carrier synchronization implementation method.

The inventive method may further comprise the steps:

Step 1, the baseband signal r of removal after down-conversion _kIn modulation intelligence, obtain one only by single-carrier signal h that additive noise polluted _k, concrete grammar is following:

A, for burst communication, utilize frame head to estimate, have

$h_k=r_k\exp (-j{\mathrm{\&theta;}}_k)=\exp (\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;}{f}_{d}k{T}_s+j{\mathrm{\&phi;}}_0)+{\hat{z}}_k,k=\mathrm{0,1},...L_0-1$

In the following formula, L ₀Length for frame head; θ _kBe the phase place of constellation point symbol of transmission, its value is θ _k∈ 0,2 π/M ..., 2 (M-1) π/M}, M is an order of modulation; f _dIt is residual frequency departure; T _sThe is-symbol cycle; φ ₀It is initial skew; J is an imaginary unit; K representes the label of k sampled point; ${\hat{z}}_k=z_k\mathrm{Exp}(-j{\mathrm{\&theta;}}_k)$ Be the additivity white complex gaussian noise, average is 0, and variance is σ ²

B, for continuous communiction, earlier with r _kBe rewritten as

r _k＝A _kexp[j(θ _k+2πf _dkT _s+φ ₀+α _k)]

In the following formula, A _kBe the instantaneous amplitude that obtains owing to noise effect, α _kIt is owing to noise effect and additional phase noise; Then, have

h _k＝|r _k|exp[jM?arg(r _k)]

Step 2, adopt searching algorithm, obtain the single-carrier signal h that obtains through step 1 based on dichotomy _kFrequency f _dWith first phase φ ₀Maximum likelihood estimator

With

Step 3, the maximum likelihood estimator that utilizes step 2 to obtain With

To baseband signal r _kCompensate, accomplish carrier synchronization and obtain demodulation result d _k

In the process of demodulation, the repeated using step 1 is to step 3, to frequency deviation f _dWith skew φ _kCarry out Tracking Estimation, thereby realize the carrier synchronization tracking of whole signals transmission.

Beneficial effect

The inventive method contrasts prior art through adopting frequency deviation and deviation estimation algorithm mutually, not only have good estimation performance, and operand is also very little, and very wide frequency offset estimation range is arranged.

Description of drawings

Fig. 1 is a burst communication data format sketch map;

Fig. 2 is a continuous communiction data format sketch map;

Fig. 3 is the schematic diagram of section's Stas ring;

Fig. 4 is the basic realization schematic diagram based on frequency offset estimating and compensation carrier synchronization of prior art;

Fig. 5 is the mould value exemplary plot of X (f);

Fig. 6 is the search procedure sketch map of frequency offset estimating value;

Fig. 7 is the performance sketch map of normalization frequency offset estimating variance

;

Fig. 8 is the performance sketch map of skew estimate variance

.

Embodiment

Below in conjunction with accompanying drawing and embodiment the inventive method is described further.

A kind of MKSP system carrier method for synchronous based on maximal possibility estimation that the present invention proposes may further comprise the steps:

Step 1, the baseband signal r of removal after down-conversion _kIn modulation intelligence, obtain a single-carrier signal h who is only polluted by additive noise _kMethod is following:

A, for burst communication, utilize frame head to estimate, because frame head data is known for receiving terminal, so have

$h_k=r_k\exp (-j{\mathrm{\&theta;}}_k)=\exp (\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;}{f}_{d}k{T}_s+j{\mathrm{\&phi;}}_0)+{\hat{z}}_k,k=\mathrm{0,1},...L_0-1---\left(8\right)$

In the formula (8), L ₀Length for frame head; θ _kBe the phase place of constellation point symbol of transmission, its value is θ _k∈ 0,2 π/M ..., 2 (M-1) π/M}, M is an order of modulation; f _dIt is residual frequency departure; T _sThe is-symbol cycle; φ ₀It is initial skew; J is an imaginary unit; K representes the label of k sampled point; ${\hat{z}}_k=z_k\mathrm{Exp}(-j{\mathrm{\&theta;}}_k)$ Be the additivity white complex gaussian noise, average is 0, and variance is σ ²

B, for continuous communiction, the data of transmission are unknown for receiving terminal, if will remove modulation intelligence, at first with r _kBe rewritten as

r _k＝A _kexp[j(θ _k+2πf _dkT _s+φ ₀+α _k)] (9)

In the formula (9), A _kBe the instantaneous amplitude that obtains owing to noise effect, α _kIt is owing to noise effect and additional phase noise.Afterwards, have

$h_k=\left|r_k\right|\exp \left[\mathrm{jM}\arg \left(r_k\right)\right]$

$=A_k\exp \left(\mathrm{jM}{\mathrm{\&theta;}}_k\right)\exp \left[j(2\mathrm{\&pi;}{\mathrm{Mf}}_{d}k{T}_s+M{\mathrm{\&phi;}}_0+M{\mathrm{\&alpha;}}_k)\right]---\left(10\right)$

$=\exp (\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;}{\mathrm{Mf}}_{d}k{T}_s+\mathrm{jM}{\mathrm{\&phi;}}_0)+{\hat{z}}_k$

In the formula (10),

is equivalent noise.It is thus clear that h _kIn do not comprised modulation intelligence.

Step 2, adopt searching algorithm, obtain the single-carrier signal h that obtains through step 1 based on dichotomy _kFrequency f _dWith first phase φ ₀Maximum likelihood estimator

With

At first, when noise is additive white Gaussian noise, obtain single-carrier signal h _kFrequency f _dWith first phase φ ₀The maximal possibility estimation expression formula.

Can find out from formula (8) and formula (10), for burst communication and continuous communiction, h _kExpression formula and parameter be that similarly all by the form of the single carrier of noise pollution, the parameter that estimate is a frequency f _dWith first phase φ ₀With the burst communication is example, provides f in the formula (8) _dAnd φ ₀The maximal possibility estimation expression formula:

Because

Be that average is 0, variance is σ ²White complex gaussian noise, so the likelihood function of formula (8) data model does

$L({\mathrm{\&phi;}}_0,f_d)\&Proportional;p(h;{\mathrm{\&phi;}}_0,f_d)$

$=\frac{1}{{\mathrm{\&pi;}}^{{\mathrm{<figure-callout\; id="0"\; label="L"\; filenames="HA20178328200910087839201E00021.png,HA20178328200910087839201E00022.png"\; state="\{\{state\}\}">L</figure-callout>}}_0}{\mathrm{\&sigma;}}^{2{\mathrm{<figure-callout\; id="0"\; label="L"\; filenames="HA20178328200910087839201E00021.png,HA20178328200910087839201E00022.png"\; state="\{\{state\}\}">L</figure-callout>}}_0}}\exp [-\frac{1}{{\mathrm{\&sigma;}}^2}{(h-x)}^H(h-x)]---\left(11\right)$

In the formula (11), $h={[{\mathrm{<figure-callout\; id="0"\; label="h"\; filenames="HA20178328200910087839201E00021.png,HA20178328200910087839201E00022.png"\; state="\{\{state\}\}">h</figure-callout>}}_0,{\mathrm{<figure-callout\; id="1"\; label="h"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">h</figure-callout>}}_1,...,h_{L_0-1}]}^T$ For receiving data vector; $x={[x_0,x_1,...,x_{L_0-1}]}^T$ Be multiple sinusoidal signal vector, wherein x _k=exp (j2 π f _dKT _s+ j φ ₀); P (h; φ ₀, f _d) be the probability density function of vectorial h, parameter wherein is φ ₀And f _dφ ₀And f _dMaximum likelihood estimator

With

Should make likelihood function L (φ ₀, f _d) value maximum, i.e. (h-x) ^H(h-x) value is minimum, can get φ in view of the above ₀And f _dThe maximal possibility estimation expression formula do

${\hat{f}}_d=\underset{f}{\arg \max }\left\{\right|\underset{k=0}{\overset{L_0-1}{\mathrm{\&Sigma;}}}{h}_k\exp (-\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;fk}{T}_s)\left|\right\}---\left(12\right)$

${\widehat{\mathrm{\&phi;}}}_0=\arg \left\{\underset{k=0}{\overset{L_0-1}{\mathrm{\&Sigma;}}}{h}_k\exp (-\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;}{\hat{f}}_{d}k{T}_s)\right\}---\left(13\right)$

Analyze in the face of above-mentioned conclusion down.Can find out from formula (13); is based on

, so its estimated performance will receive the influence of frequency offset estimating error.Utilizing formula (7) to accomplish carrier synchronization, realize in the process of demodulation, because

Inaccurate and phase error accumulation also can be increasing, and then influence demodulation result, so need be constantly to f in the process of demodulation _dAnd φ _kEstimate, promptly frequency deviation followed the tracks of that the method for estimation when estimation approach is with continuous communiction during tracking is identical.When the influencing of consideration of noise not, formula (13) can be rewritten as

${\widehat{\mathrm{\&phi;}}}_0=\arg \left\{\underset{k=0}{\overset{L_0-1}{\mathrm{\&Sigma;}}}\exp (\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;}{f}_{d}k{T}_s+j{\mathrm{\&phi;}}_0)\exp (-\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;}{\hat{f}}_{d}k{T}_s)\right\}$

$=\arg \left\{\underset{k=0}{\overset{L_0-1}{\mathrm{\&Sigma;}}}\exp \left(j{\mathrm{\&phi;}}_0\right)\exp \right[\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;}(f_d-{\hat{f}}_d)k{T}_s]\}---\left(14\right)$

$=\arg \left\{\exp \right[j{\mathrm{\&phi;}}_0+\mathrm{j\&pi;}(f_d-{\hat{f}}_d)(L_0-1)T_s]$

$\&times;\underset{k=0}{\overset{L_0-1}{\mathrm{\&Sigma;}}}\exp \left[\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;}(f_d-{\hat{f}}_d)(k-(L_0-1)/2)T_s\right]\}$

Obviously, and likes $\underset{k=0}{\overset{L_0-1}{\mathrm{\&Sigma;}}}\mathrm{Exp}\left[\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;}(f_d-{\hat{f}}_d)(k-(L_0-1)/2)T_s\right]$ Be an arithmetic number, so formula (14) can be abbreviated as

${\widehat{\mathrm{\&phi;}}}_0={\mathrm{\&phi;}}_0+\mathrm{\&pi;}(f_d-{\hat{f}}_d)(L_0-1)T_s---\left(15\right)$

When

Depart from f _dWhen too many, the phase estimation error will become greatly, and

Not not have partially to estimate.Convolution (15), the first phase estimated value of frame head intermediate symbols does

${\widehat{\mathrm{\&phi;}}}_{(L_0-1)/2}={\widehat{\mathrm{\&phi;}}}_0+\mathrm{\&pi;}{\hat{f}}_d(L_0-1)T_s={\mathrm{\&phi;}}_0+\mathrm{\&pi;}{f}_d(L_0-1)T_s={\mathrm{\&phi;}}_{(L_0-1)/2}---\left(16\right)$

,

estimate partially that the variance that therefore can use

is as the standard of weighing the skew estimated performance it is thus clear that being the nothing of .

For continuous communiction, ask φ ₀And f _dThe method of maximum likelihood estimator identical during with burst communication, the maximal possibility estimation expression formula does

${\hat{f}}_d=\frac{1}{M}\underset{f}{\arg \max }\left\{\right|\underset{k=0}{\overset{K-1}{\mathrm{\&Sigma;}}}{h}_k\exp (-\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;fk}{T}_s)\left|\right\}---\left(17\right)$

${\widehat{\mathrm{\&phi;}}}_0^{\&prime;}=\frac{1}{M}\arg \left\{\underset{k=0}{\overset{K-1}{\mathrm{\&Sigma;}}}{h}_k\exp (-\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;M}{\hat{f}}_{d}k{T}_s)\right\}---\left(18\right)$

In formula (17), the formula (18), K is the data length that is used to estimate, it is relevant with needed estimation accuracy.Owing to be mpsk signal, estimation is come out

M phase ambiguity arranged, that is, actual skew estimated value should be ${\widehat{\mathrm{\&phi;}}}_0^{\&prime;}+2\mathrm{\&pi;\; m}/M,m=\mathrm{0,1},...,M-1$ In one, for any algorithm, all can't remove the phase ambiguity of blind mpsk signal, can only remove phase ambiguity through known unique word.The demodulation result at unique word place when at first not removed phase ambiguity $d_k=r_k\mathrm{Exp}(-\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;}{\hat{f}}_{d}k{T}_s-j{\widehat{\mathrm{\&phi;}}}_0^{\&prime;}),k=\mathrm{0,1},...L_s-1,$ The data content of unique word is s _k, k=0,1 ... L _s-1, then can remove phase ambiguity by through type (19), obtain actual first phase maximal possibility estimation expression formula:

${\widehat{\mathrm{\&phi;}}}_0=\underset{\mathrm{\&phi;}}{\arg \min }\{\mathrm{\&phi;}-\arg \left(\underset{k=0}{\overset{L_s-1}{\mathrm{\&Sigma;}}}{d}_k{s}_k^{*}\right)\},\mathrm{\&phi;}\&Element;\{{\widehat{\mathrm{\&phi;}}}_0^{\&prime;}+2\mathrm{\&pi;m}/M,m=\mathrm{0,1},...,M-1\}---\left(19\right)$

It is the same with the situation of burst communication,

Also be φ ₀Inclined to one side estimation arranged, but through type (16) obtains

Then be

Nothing estimate partially, therefore can be used as the standard of weighing the skew estimated performance.

At this moment, utilize searching algorithm, obtain f fast based on dichotomy _dAnd φ ₀Maximum likelihood estimator

With

Process is following:

Can find out f from formula (12) and formula (17) _dMaximum likelihood estimator

Analytic solutions can not directly not calculate, and therefore need to adopt the fast search algorithm based on dichotomy.Because for burst communication and continuous communiction; The expression formula of

is similarly, is the statement that example is carried out this searching algorithm with the burst communication therefore.

Definition

$X\left(f\right)=\underset{k=0}{\overset{L_0-1}{\mathrm{\&Sigma;}}}{h}_k\exp (-\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;fk}{T}_s)---\left(20\right)$

Obviously, f _dMaximum likelihood estimator Be the maximum frequency of mould value that makes X (f).When consideration of noise

not, the mould value of X (f) can be written as

$\left|X\left(f\right)\right|=\left|\underset{k=0}{\overset{L_0-1}{\mathrm{\&Sigma;}}}{x}_k\exp (-\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;fk}{T}_s)\right|=\left|\frac{\sin \left[\mathrm{\&pi;}(f-f_d)L_0{T}_s\right]}{\sin \left[\mathrm{\&pi;}(f-f_d)T_s\right]}\right|---\left(21\right)$

As (f-f _d) T _s, sin [π (f-f is arranged at＜＜1 o'clock _d) T _s] ≈ π (f-f _d) T _s, then formula (21) can be similar to and be written as

|X(f)|＝L ₀|sinc[π(f-f _d)L ₀T _s]| (22)

Sinc (x) is defined as sin (x)/x in the formula (22).The mould value example of X (f) is as shown in Figure 5.

At first confirm searching times, method is following:

In reality, f _dCertain scope is all arranged, promptly | f _d|＜f _Max, f _MaxIt is possible maximum frequency deviation.Needed searching times and f _MaxRelevant with desired Frequency Estimation precision.

Work as f _Max≤1/L ₀T _sThe time, if the maximum search number of times is N, the error between the estimated value that then obtains according to this searching algorithm and the maximum likelihood estimator of reality is less than f _Max/ 2 ^NIn the digital receiver of practical communication system, normalized frequency fT _sAll represent, if fT with fixed-point number _sBe quantized into the n bit, the normalized frequency of the minimum that then can represent is 1/2 ⁿ, this moment, needed maximum search number of times did

$N=\mathrm{ceil}\left({\log }_2^{f_{\max }{T}_s\&times;2^n}\right)=n+\mathrm{ceil}\left({\log }_2^{f_{\max }{T}_s}\right)---\left(23\right)$

Thus it is clear that, maximum search times N and maximum frequency deviation f _MaxBecome logarithmic relationship, explain that the operand of this searching algorithm is very little.

Work as f _Max＞1/L ₀T _sThe time, need when initial ranging, need calculate 2t X (f) through an initial ranging, through after the initial ranging, be equivalent to the maximum frequency deviation scope is reduced into f _Max=1/2L ₀T _s, in the subsequent searches process, needed maximum search number of times does

$N=n+\mathrm{ceil}\left({\log }_2^{f_{\max }{T}_s}\right)=n-\mathrm{ceil}\left({\log }_2^{L_0}\right)-1---\left(24\right)$

After confirming searching times, adopt searching algorithm to obtain f based on dichotomy _dAnd φ ₀Maximum likelihood estimator

With

As can be seen from Figure 5, | X (f) | the main value width be 2/L ₀T _s, therefore, suppose f earlier _Max≤1/L ₀T _s, then the step of searching algorithm is:

Frequency values f in the middle of I, two search rates in the time of will at every turn searching for _MidFrequency step Δ f during with each search is initialized as

f _mid＝0，Δf＝f _max/2 (25)

II, make f ₁=f _Mid-Δ f and f ₂=f _Mid+ Δ f calculates X (f according to formula (20) then ₁) and X (f ₂), relatively | X (f ₁) | with | X (f ₂) |, if | X (f ₁) |＞| X (f ₂) |, then make f _Mid=f ₁Otherwise make f _Mid=f ₂If reached the maximum search number of times, then jump to Step II I; Otherwise make Δ f=Δ f/2, repeating step II.

III, the f after will searching for for the last time _MidEstimated value as frequency deviation

Through type (13) calculates then

Obtain by formula (16) again

Work as f _Max＞1/L ₀T _sThe time, at the beginning need be through an initial ranging.Schilling t=ceil (f _MaxL ₀T _s), function ceil (x) expression rounds up, and makes f then _i=(i-1)/L ₀T _s+ 1/2L ₀T _s, i=-t+1 ,-t+2 ..., t, and calculate X (f according to formula (20) _i), relatively more all | X (f _i) |, will make | X (f _i) | maximum frequency is designated as f _{I max}, make f at last _Mid=f _{I max}With Δ f=1/4L ₀T _s, forward Step II to and continue search procedure.

The search procedure of frequency offset estimating value is as shown in Figure 6, and search procedure has been shown among the figure 3 times.

If continuous communiction then only needs when Step II I order ${\hat{f}}_d=f_{\mathrm{Mid}}/M,$ Through type (19) calculates then

Step 3, the maximum likelihood estimator that utilizes step 2 to obtain With To baseband signal r _kCompensate, accomplish carrier synchronization and obtain demodulation result d _k:

$d_k=r_k\exp (-\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">j</figure-callout>}2\mathrm{\&pi;}{\hat{f}}_{d}k{T}_s-j{\widehat{\mathrm{\&phi;}}}_0)---\left(26\right)$

In the process of demodulation, the repeated using step 1 is to step 3, to f _dAnd φ _kCarry out Tracking Estimation, thereby realize the carrier synchronization tracking of whole signals transmission.

Weighing the standard of an algorithm for estimating performance quality, promptly is the size of seeing its variance.Because the used algorithm for estimating of the present invention is maximal possibility estimation in essence, therefore the variance of resulting estimator promptly is a carat Mei-Luo (Cramer-Rao) boundary.Carat Mei-Luo circle is the lower limit of all unbiased estimator variances, and therefore algorithm for estimating of the present invention has optimum performance.

The normalized frequency offset estimation variance

and the phase offset estimation variance

as a measure of the frequency offset estimate

and the phase offset estimator

performance criteria.Carat Mei-Luo circle of frequency offset estimating amount

does

$\mathrm{CRB}\left({\hat{f}}_d\right)=\frac{1}{-E\left[\frac{{\&PartialD;}^2\ln p(h;f_d)}{{\&PartialD;f}_d^2}\right]}=\frac{6}{{\left(2\mathrm{\&pi;}\right)}^2{L}_0({{\mathrm{<figure-callout\; id="0"\; label="L"\; filenames="HA20178328200910087839201E00021.png,HA20178328200910087839201E00022.png"\; state="\{\{state\}\}">L</figure-callout>}}_0}^2-1)\left(\mathrm{SNR}\right)}---\left(27\right)$

In the formula (27) $\mathrm{SNR}=\frac{1}{{\mathrm{\&sigma;}}^2}$ Be signal to noise ratio (signal power is normalized to 1).Carat Mei-Luo circle of phase biased estimator

does

$\mathrm{CRB}\left({\widehat{\mathrm{\&phi;}}}_{(L_0-1)/2}\right)=\frac{1}{-e\left[\frac{{\&PartialD;}^2\ln p(h;{\mathrm{\&phi;}}_{(L_0-1)/2})}{\&PartialD;{\widehat{\mathrm{\&phi;}}}_{(L_0-1)/2}^2}\right]}=\frac{1}{{2L}_0\left(\mathrm{SNR}\right)}---\left(28\right)$

But the present invention removes to approach maximal possibility estimation with searching algorithm, so the performance of frequency offset estimating also has relation with the number of times of searching for.If hypothesis f _Max=1/ (L ₀T _s), the maximum search number of times is N, then frequency resolution (minimum frequency that can represent) is 1/ (2 ^NL ₀T _s), therefore have

So $E\left[{\left(({\hat{f}}_d-f_d)T_s\right)}^2\right]\&GreaterEqual;1/{\left(2^N{L}_0\right)}^2=1/\left(4^{\mathrm{<figure-callout\; id="0"\; label="N\; L"\; filenames="HA20178328200910087839201E00021.png,HA20178328200910087839201E00022.png"\; state="\{\{state\}\}">N</figure-callout>}}{L}_{}^{}{}_0^2\right).$ Therefore, the normalization frequency offset estimating variance of searching algorithm of the present invention

Lower limit satisfy

$E\left[{\left(({\hat{f}}_d-f_d)T_s\right)}^2\right]\&GreaterEqual;\max \{\mathrm{CRB}\left({\hat{f}}_d\right),1/\left(4^{\mathrm{<figure-callout\; id="0"\; label="N\; L"\; filenames="HA20178328200910087839201E00021.png,HA20178328200910087839201E00022.png"\; state="\{\{state\}\}">N</figure-callout>}}{L}_{}^{}{}_0^2\right)\}---\left(29\right)$

In order to check the performance of above-mentioned algorithm, the present invention has done the lots of emulation test, emulation at different length L ₀The performance of following algorithm of the present invention, and compare with carat Mei-Luo lower bound.Monte Carlo (Monte-Carlo) method is adopted in emulation, and simulation times is 5000 times, and the parameter in the emulation is f _Max=1/ (L ₀T _s), the maximum search number of times is N=5.

Provided normalization frequency offset estimating variance among Fig. 7 Performance, as can be seen from the figure, when $\mathrm{CRB}\left({\hat{f}}_d\right)>1/\left(4^{\mathrm{<figure-callout\; id="0"\; label="N\; L"\; filenames="HA20178328200910087839201E00021.png,HA20178328200910087839201E00022.png"\; state="\{\{state\}\}">N</figure-callout>}}{L}_{}^{}{}_0^2\right)$ The time, the very approaching carat Mei-Luo lower bound of the performance of searching algorithm of the present invention, otherwise near 1/ (4 ^NL ₀ ²), this is consistent with theory analysis.

Provided the performance of skew estimate variance

among Fig. 8; As can be seen from the figure, the performance of algorithm of the present invention is very near a carat Mei-Luo lower bound.The low and L when signal to noise ratio ₀More in short-term, the performance Bick Latin America-Luo lower bound of algorithm of the present invention has deterioration slightly, and this is that the frequency offset estimating error becomes bigger, thereby has worsened the skew estimation performance because in this case.Under the low signal-to-noise ratio condition, this can be through lengthening L ₀Solve.

Embodiment

Suppose certain BPSK modulating system, character rate is f _s=1/T _s=1MBaud, maximum frequency deviation is f _Max=20kHz.For burst mode, frame head is long to be L ₀=128; For continuous mode, suppose to be used to estimate that frequency deviation and the used data length of skew also are K=L ₀=128.Then the practical implementation step of carrier synchronization is following:

Step 1, the baseband signal r of removal after down-conversion _kIn modulation intelligence, obtain a single-carrier signal h who is only polluted by additive noise _k

For the burst communication pattern, frame head is known for receiving terminal, establishes the local frame head of receiving terminal and is s through BPSK mapping back _k, k=0,1 ..., L0-1, the header signal of then removing behind the modulation intelligence does

$h_k=r_k{s}_k^{*},k=\mathrm{0,1},...,L_0-1---\left(30\right)$

For the continuous communiction pattern, data are unknown for receiving terminal, and then receiving data can remove modulation intelligence by through type (31)

h _k＝|r _k|exp[j2arg(r _k)] (31)

Step 2, adopt searching algorithm, obtain the single-carrier signal h that obtains through step 1 based on dichotomy _kFrequency f _dWith first phase φ ₀Maximum likelihood estimator

With

Obtaining φ ₀And f _dThe maximal possibility estimation expression formula after, confirm earlier the maximum search number of times, because 1/L ₀T _s=1000kHz/128=7.8125kHz＜20kHz=f _MaxSo, need carry out an initial ranging, because t=ceil (f _MaxL ₀T _s)=ceil (20*128/1000)=3 is so need to calculate f _i=(i-1)/L ₀T _s+ 1/2L ₀T _s=(i-1/2) * 7.8125kHz, i=2 ,-1 ..., the X (f at 3 these 6 frequency places _i), and relatively more all | X (f _i) |, made | X (f _i) | maximum frequency f _{I max}If supposition is when realizing, normalized frequency fT _sRepresent with 14 bit fixed point numbers, then calculate and to get the maximum search number of times by formula (24) $N=14-\mathrm{Ceil}\left({\mathrm{<figure-callout\; id="2"\; label="Log"\; filenames="HA20178328200910087839201E00012.png,HA20178328200910087839201E00021.png"\; state="\{\{state\}\}">Log</figure-callout>}}_2^{128}\right)-1=6.$ Search step is following:

I, at first make f _Mid=f _Imax, Δ f=1/4L ₀T _s=1.953125kHz jumps to Step II.

II, make f ₁=f _Mid-Δ f and f ₂=f _Mid+ Δ f calculates X (f according to formula (20) then ₁) and X (f ₂), relatively | X (f ₁) | with | X (f ₂) |, if | X (f ₁) |＞| X (f ₂) |, then make f _Mid=f ₁Otherwise make f _Mid=f ₂If reached maximum search number of times 6, then jump to Step II I; Otherwise make Δ f=Δ f/2, repeating step II.

III, for burst communication, with the 6th time the search after f _MidEstimated value as frequency deviation

Through type (13) calculates then

The process of accomplishing search and estimating; And for continuous communiction, then have ${\hat{f}}_d=f_{\mathrm{Mid}}/2,$ Through type (19) obtains then

The process of accomplishing search and estimating.

Step 3, the maximum likelihood estimator that utilizes step 2 to obtain

With

To baseband signal r _kCompensate, accomplish carrier synchronization and obtain demodulation result d _k:

Utilize formula (26) to accomplish carrier synchronization and demodulation, in the process of demodulation, the repeated using step 1 is to step 3, to f _dAnd φ _kCarry out Tracking Estimation, thereby realize the carrier synchronization tracking of whole signals transmission.

Claims (1)

1. MPSK system carrier method for synchronous based on maximal possibility estimation is characterized in that may further comprise the steps:

Step 1, the baseband signal r of removal after down-conversion _kIn modulation intelligence, obtain one only by single-carrier signal h that additive noise polluted _k, concrete grammar is following:

A, for burst communication, utilize frame head to estimate, have

$h_k=r_k\exp (-j{\mathrm{\&theta;}}_k)=\exp (j2\mathrm{\&pi;}{f}_{d}k{T}_s+j{\mathrm{\&phi;}}_0)+{\hat{z}}_k,$ k＝0，1，…L ₀-1

In the following formula, L ₀Length for frame head; θ _kBe the phase place of constellation point symbol of transmission, its value is θ _k∈ 0,2 π/M ..., 2 (M-1) π/M}, M is an order of modulation; f _dIt is residual frequency departure; T _sThe is-symbol cycle; φ ₀It is initial skew; J is an imaginary unit; K representes the label of k sampled point;

Be the additivity white complex gaussian noise, average is 0, and variance is σ ²

B, for continuous communiction, earlier with r _kBe rewritten as

r _k＝Akexp[j(θ _k+2πf _dkT _s+φ ₀+α _k)]

In the following formula, A _kBe the instantaneous amplitude that obtains owing to noise effect, α _kIt is owing to noise effect and additional phase noise; Then, have

h _k＝|r _k|exp[jMarg(r _k)]

Step 2, adopt searching algorithm, obtain the single-carrier signal h that obtains through step 1 based on dichotomy _kFrequency f _dWith first phase φ ₀Maximum likelihood estimator

With

The method that this step realizes is following:

At first, when noise is additive white Gaussian noise, obtain single-carrier signal h _kFrequency f _dWith first phase φ ₀The maximal possibility estimation expression formula:

(1) for burst communication, frequency f _dWith first phase φ ₀The maximal possibility estimation expression formula be respectively

${\hat{f}}_d=\underset{f}{\arg \max }\left\{\right|\underset{k=0}{\overset{L_0-1}{\mathrm{\&Sigma;}}}{h}_k\exp (-j2\mathrm{\&pi;fk}{T}_s)\left|\right\}$

${\widehat{\mathrm{\&phi;}}}_0=\arg \left\{\underset{k=0}{\overset{L_0-1}{\mathrm{\&Sigma;}}}{h}_k\exp (-j2\mathrm{\&pi;}{\hat{f}}_{d}k{T}_s)\right\}$

(2) for continuous communiction, frequency f _dWith first phase φ ₀The maximal possibility estimation expression formula be respectively

${\hat{f}}_d=\underset{f}{\frac{1}{M}\arg \max }\left\{\right|\underset{k=0}{\overset{K-1}{\mathrm{\&Sigma;}}}{h}_k\exp (-j2\mathrm{\&pi;fk}{T}_s)\left|\right\}$

${\widehat{\mathrm{\&phi;}}}_0=\underset{\mathrm{\&phi;}}{\arg \min }\{\mathrm{\&phi;}-\arg \left(\underset{k=0}{\overset{L_s-1}{\mathrm{\&Sigma;}}}{d}_k{s}_k^{*}\right)\},$ $\mathrm{\&phi;}\&Element;\{{\widehat{\mathrm{\&phi;}}}_0^{\&prime;}+2\mathrm{\&pi;m}/M,m=\mathrm{0,1},...,M-1\}$

${\widehat{\mathrm{\&phi;}}}_0^{\&prime;}=\frac{1}{M}\arg \left\{\underset{k=0}{\overset{K-1}{\mathrm{\&Sigma;}}}{h}_k\exp (-j2\mathrm{\&pi;M}{\hat{f}}_{d}k{T}_s)\right\}$

Wherein, K is the data length that is used to estimate,

Be the maximal possibility estimation of first phase that ambiguity is arranged, L _sBe unique word length, d _kThe demodulation result at unique word place when not removing phase ambiguity, S _kData content for unique word;

Afterwards, adopt searching algorithm to obtain f based on dichotomy _dAnd φ ₀Maximum likelihood estimator

With

Detailed process is following:

(1) at first confirms searching times

Work as f _Max≤1/L ₀T _sThe time, needed maximum search number of times does

$N=\mathrm{ceil}\left({\log }_2^{f_{\max }{T}_s\&times;2^n}\right)=n+\mathrm{ceil}\left({\log }_2^{f_{\max }{T}_s}\right)$

Work as f _Max＞1/L ₀T _sThe time, earlier through an initial ranging, the maximum frequency deviation scope is reduced into f _Max=1/2L ₀T _s, in the subsequent searches process, needed maximum search number of times does

$N=n+\mathrm{ceil}\left({\log }_2^{f_{\max }{T}_s}\right)=n-\mathrm{ceil}\left({\log }_2^{L_0}\right)-1$

Wherein, f _MaxBe possible maximum frequency deviation, the bit number of n for quantizing;

(2) confirm searching times after, adopt searching algorithm to obtain f based on dichotomy _dAnd φ ₀Maximum likelihood estimator

With

For burst communication, work as f _Max≤1/L ₀T _sThe time, the step of searching algorithm is:

Frequency values f in the middle of I, two search rates in the time of will at every turn searching for _MidFrequency step Δ f during with each search is initialized as

f _mid＝0，Δf＝f _max/2

II, make f ₁=f _Mid-Δ f and f ₂=f _Mid+ Δ f, basis then

Calculate X (f ₁) and X (f ₂), relatively | X (f ₁) | with | X (f ₂) |, if | X (f ₁) |＞| X (f ₂) |, then make f _Mid=f ₁Otherwise make f _Mid=f ₂If reached the maximum search number of times, then jump to Step II I; Otherwise make Δ f=Δ f/2, repeating step II;

III, the f after will searching for for the last time _MidEstimated value as frequency deviation

Pass through then ${\widehat{\mathrm{\&phi;}}}_0=\mathrm{Arg}\left\{\underset{k=0}{\overset{L_0-1}{\mathrm{\&Sigma;}}}{h}_k\mathrm{Exp}(-j2\mathrm{\&pi;}{\hat{f}}_{d}k{T}_s)\right\}$ Calculate Again by ${\widehat{\mathrm{\&phi;}}}_{(L_0-1)/2}={\widehat{\mathrm{\&phi;}}}_0+\mathrm{\&pi;}{\hat{f}}_d(L_0-1)T_s={\mathrm{\&phi;}}_0+\mathrm{\&pi;}{f}_d(L_0-1)T_s={\mathrm{\&phi;}}_{(L_0-1)/2}$ Obtain

Work as f _Max＞1/L ₀T _sThe time, at the beginning need be through an initial ranging, shilling f=ceil (f _MaxL ₀T _s), function ceil (x) expression rounds up, and makes f then _i=(i-1)/L ₀T _s+ 1/2L ₀T _s, i=-t+1 ,-t+2 ..., t, and according to

Calculate X (f _i), relatively more all | X (f _i) |, will make | X (f _i) | maximum frequency is designated as f _{I max}, make f at last _Mid=f _{I max}With Δ f=1/4L ₀T _s, forward Step II to and continue search procedure;

If continuous communiction then only needs when Step II I, will Pass through then ${\widehat{\mathrm{\&phi;}}}_0=\underset{\mathrm{\&phi;}}{\mathrm{Arg}\mathrm{Min}}\{\mathrm{\&phi;}-\mathrm{Arg}\left(\underset{k=0}{\overset{L_s-1}{\mathrm{\&Sigma;}}}\right)d_k{s}_k^{*}\},$ $\mathrm{\&phi;}\&Element;\{{\widehat{\mathrm{\&phi;}}}_0^{\&prime;}+2\mathrm{\&pi;\; m}/M,m=\mathrm{0,1},...,M-1\}$ Calculate

Step 3, the maximum likelihood estimator that utilizes step 2 to obtain

With

To baseband signal r _kCompensate, accomplish carrier synchronization and obtain demodulation result d _k

In the process of demodulation, the repeated using step 1 is to step 3, to frequency deviation f _dWith skew φ _kCarry out Tracking Estimation, thereby realize the carrier synchronization tracking of whole signals transmission.

Images