DSPRelated.com
Forums

Software PLL - the math that i did looks too simple to be true - please confrim

Started by mir_aculous February 7, 2011
>On 02/07/2011 09:50 PM, mir_aculous wrote: >>> On Feb 7, 8:39=A0pm, Tim Wescott<t...@seemywebsite.com> wrote: >>>> On 02/07/2011 02:44 PM, robert bristow-johnson wrote: >>>> >>>> >>>>> since the input to the NCO directly controls frequency, you need to >>>>> model (and it's simple, just an integrator) the transfer function >> from >>>>> frequency to phase, since it's (usually) the phase difference that >>>>> comes outa the phase discriminator. =A0the integrator that is >> inherent >>>>> to the NCO usually obviates the need for the "P" in the PID >>>>> controller. >>>> >>>> Eh? =A0Most phase-locked loops that I do use a PI controller. >>> >>> ooops. i meant that in the actual guts of the controller, you would >>> lose the "I", because of the inherent I in the NCO (if what comes outa >>> the phase discriminator is the difference of phase, not the difference >>> of frequency). >>> >>>> =A0One rarely >>>> needs differential, but one occasionally needs an extra low-pass
filter
>>>> if you really want to clean up jitter. >>>> >>> >>> what i meant was this identity: >>> >>> >>> sin(w0+phi)*cos(w0) =3D 0.5*sin(phi) + 0.5*sin(2*w0+phi) >>> >>> because the integrator might take care of it, you might not need the >>> LPF, but we want the sin(phi) term to survive as it is pretty close to >>> phi and measures our phase difference. >>> >>>>> another way to mathematically model a phase discriminator is to look >>>>> at it as the LPFed output of a multiplier between the two waveforms >>>>> that are virtually the same frequency. =A0when the PLL is "locked", >> tha= >>> t >>>>> means the two waveforms will be 90 degrees outa phase. >>>> >>>>> of course, when they get more and more out of phase, there is a non- >>>>> linearity in the phase discriminator, if you do it this way. >>>> >>>>> if you need your oscillator output to be in the same phase as the >>>>> input, then use the same phase going into the NCO for two quadrature >>>>> (that is 90 degrees outa phase with each other) table lookups. >>>> >>>> That depends on the phase detector. =A0If you're locking to edges >> (which >>>> is fairly easy to do with the right counter resources) then you can >> make >>>> your phase detector act pretty much like a regular old
phase/frequency
>>>> detector. =A0If you're locking to data, then you're locking to data
and
>>>> limited by the data stream -- and that's not going to be like locking >> to >>>> a sinusoid. >>> >>> agreed. but multiplying two square waves together (and LPFing) can >>> get that same phi term that goes to zero when the two are 90 degrees >>> outa phase. i seem to remember seeing an old design (i wonder if it >>> *was* a 4016 or similar) where the phase discriminator was essentially >>> an XOR gate. it's sorta similar. >>> >>> r b-j >>> >>> >>> >> >> Before, I get confused further about the implementation process, I want
to
>> tell you all that my ultimate goal is to track the carrier and the code
in
>> a direct sequence spread spectrum signal. >> >> First, I want to implement the basic PLL where I can input a signal such
as
>> an FM signal and try to lock on to the carrier. If that works then I
can
>> move on towards my goal of locking on to the carrier of a DSSS signal
and
>> then doing clock/code tracking by implementing a code tracking loop as >> given in many gps related books. But first, the basic loop to track a >> sinusoid. >> >>> since the input to the NCO directly controls frequency, you need to >>> model (and it's simple, just an integrator) the transfer function from >>> frequency to phase, since it's (usually) the phase difference that >>> comes outa the phase discriminator. the integrator that is inherent >>> to the NCO usually obviates the need for the "P" in the PID >>> controller. >> >> r-b-j, can you explain a little more about the direct implementation
that
>> you are talking about. This is what I understand so far. >> >> phase_error = LPF ( input sample * VCO sample ), which makes sense as
it
>> ultimately splits down to a phase term and another high frequency term >> which is low pass filtered. So, once I have this phase term what do I
do
>> next? > >Robert is giving you a very limited view of phase detectors. Some do, >indeed, multiply the incoming signals (either as sinusoids, or as square >waves, in which case you find that an XOR gate acts the same as a >multiplication). Some compare edges (search on "three state phase >detector" or "phase-frequency detector"). Some PLL's are extracting >clock from data, and there are various ways to do that (search on "data >driven PLL" or "Costas loop"). In the case of DSSS, you run at least >three correlators, an 'early' one, an 'on time' one, and a 'late' one. >You compare the relative magnitudes of the outputs of the early and late >correlators and use that to correct the loop. > >> Can you put it down in a mathematical equation as it's easier to
understand
>> that way. > >You describe the behavior in terms of a difference equation, then you >take the z transform of the difference equation. If, for instance, you >have a phase-increment NCO, then at each iteration of your clock you >find that your NCO phase goes > >p_O(n) = p_O(n-1) + p_i(n) > >where p_O is the oscillator phase, and p_i is the commanded phase from >the loop filter. Taking the z transform of this you get > >P_O = P_O * z^-1 + P_i > >or > >z P_O = P_O + z P_i > >or > >P_O z >--- = ----- >P_i z - 1 >
Extending, on this and using a loop filter with transfer function shown below (assuming all the gains are set to 1 for cleaner representation), H(s) = (tau2.s+1)/(tau1.s) , I get the difference equation as u(n) = e(n)[2.tau2/Ts + 1](Ts/(2.tau1)) - e(n-1)[1 - 2.tau2/Ts](Ts/(2.tau1))+ u(n-1) where, u(n) is the loop filter output at time t=n, e(n) is the output of the phase detector, tau1, tau2 are the loop filer constants, Ts is the sampling width or 1/Fs where, Fs is the sampling frequency Finally, as you mentioned earlier u(n) is the input to the VCO thus, VCO(n) = VCO(n-1) + u(n) For DSSS, I found out that a common phase detection scheme is e(n) = arctan(Q/I) where, Q and I are outputs of the inphase and quadrature components of the prompt correlator.
>You model the phase detector by a similar process (and find that it's >just a gain, or a system that can be linearized to get a gain when the >loop is locked), then you make a loop filter that will get you your >desired loop performance. > >See http://www.wescottdesign.com/articles/zTransform/z-transforms.html >and http://www.wescottdesign.com/actfes/actfes.html. > >-- > >Tim Wescott >Wescott Design Services >http://www.wescottdesign.com > >Do you need to implement control loops in software? >"Applied Control Theory for Embedded Systems" was written for you. >See details at http://www.wescottdesign.com/actfes/actfes.html >
On Feb 8, 12:50&#4294967295;am, "mir_aculous" <gnu.fanz@n_o_s_p_a_m.gmail.com>
wrote:
> >On Feb 7, 8:39=A0pm, Tim Wescott <t...@seemywebsite.com> wrote: > >> On 02/07/2011 02:44 PM, robert bristow-johnson wrote: > > >> > since the input to the NCO directly controls frequency, you need to > >> > model (and it's simple, just an integrator) the transfer function > from > >> > frequency to phase, since it's (usually) the phase difference that > >> > comes outa the phase discriminator. =A0the integrator that is > inherent > >> > to the NCO usually obviates the need for the "P" in the PID > >> > controller. > > >> Eh? =A0Most phase-locked loops that I do use a PI controller. > > >ooops. &#4294967295;i meant that in the actual guts of the controller, you would > >lose the "I", because of the inherent I in the NCO (if what comes outa > >the phase discriminator is the difference of phase, not the difference > >of frequency). > > >> =A0One rarely > >> needs differential, but one occasionally needs an extra low-pass filter > >> if you really want to clean up jitter. > > >what i meant was this identity: > > > &#4294967295; &#4294967295;sin(w0+phi)*cos(w0) =3D 0.5*sin(phi) + 0.5*sin(2*w0+phi) > > >because the integrator might take care of it, you might not need the > >LPF, but we want the sin(phi) term to survive as it is pretty close to > >phi and measures our phase difference. > > >> > another way to mathematically model a phase discriminator is to look > >> > at it as the LPFed output of a multiplier between the two waveforms > >> > that are virtually the same frequency. =A0when the PLL is "locked", > tha= > >t > >> > means the two waveforms will be 90 degrees outa phase. > > >> > of course, when they get more and more out of phase, there is a non- > >> > linearity in the phase discriminator, if you do it this way. > > >> > if you need your oscillator output to be in the same phase as the > >> > input, then use the same phase going into the NCO for two quadrature > >> > (that is 90 degrees outa phase with each other) table lookups. > > >> That depends on the phase detector. =A0If you're locking to edges > (which > >> is fairly easy to do with the right counter resources) then you can > make > >> your phase detector act pretty much like a regular old phase/frequency > >> detector. =A0If you're locking to data, then you're locking to data and > >> limited by the data stream -- and that's not going to be like locking > to > >> a sinusoid. > > >agreed. &#4294967295;but multiplying two square waves together (and LPFing) can > >get that same phi term that goes to zero when the two are 90 degrees > >outa phase. &#4294967295;i seem to remember seeing an old design (i wonder if it > >*was* a 4016 or similar) where the phase discriminator was essentially > >an XOR gate. &#4294967295;it's sorta similar. > > >r b-j > > Before, I get confused further about the implementation process, I want to > tell you all that my ultimate goal is to track the carrier and the code in > a direct sequence spread spectrum signal. > > First, I want to implement the basic PLL where I can input a signal such as > an FM signal and try to lock on to the carrier. If that works then I can > move on towards my goal of locking on to the carrier of a DSSS signal and > then doing clock/code tracking by implementing a code tracking loop as > given in many gps related books. But first, the basic loop to track a > sinusoid. > > > since the input to the NCO directly controls frequency, you need to > > model (and it's simple, just an integrator) the transfer function from > > frequency to phase, since it's (usually) the phase difference that > > comes outa the phase discriminator. &#4294967295;the integrator that is inherent > > to the NCO usually obviates the need for the "P" in the PID > > controller. > > r-b-j, can you explain a little more about the direct implementation that > you are talking about. This is what I understand so far. > > &#4294967295; phase_error = LPF ( input sample * VCO sample ),
shall i make the assumption that both are sinusoidal? the NCO can have a sinusoidal lookup table. you can set up math for the square wave case or if one is square and the other a ramp or even if one is sinusoidal and the other is square or ramp. it really doesn't matter. if the waveforms are symmetrical and zero-mean and you define them to be "in phase" if their zero crossings line up, then this PLL will be "locked" when they're 90 out of phase. this is easy to model if both waveforms are so close in (instantaneous) frequency (the rate of change of phase) that we model them as identical in frequency with a possibly time-varying change in phase (of one relative to the other) that will account for a tiny frequency difference. the NCO output: y_nco[n] = cos( phi_nco[n] ) where phi_nco[n] = phi_nco[n-1] + G_nco*x_nco[n] and x_nco[n] is the control "voltage" input for the NCO. radians are actually dimensionless, but G_nco is the "NCO gain" and might be in units of radian/nco_ref or G_nco = 1/nco_ref. when x_nco[n] = nco_ref, the phase of the NCO advances 1 radian per sample. but if you consider x_nco[n] dimensionless, so is G_nco. the input is x[n] = A * sin( phi[n] ) where phi[n] = phi[n-1] + omega[n] at times that the frequency omega[n] is mostly constant, we can say that and, for this simple analysis, omega[n] gotta be close to G_nco*x_nco[n]. so the output of the simple 4-quadrant multiplier is x[n]*y_nco[n] = A*sin(phi[n]) * cos(phi_nco[n]) = A*sin(phi[n]) * cos(phi_nco[n]) = A/2*( sin((phi[n]+phi_nco[n])/2) + sin((phi[n]-phi_nco[n])/2) ) the first term is sin((phi[n]+phi_nco[n])/2) = sin( (omega[n]+G_nco[n]*x_nco[n])/2 + (phi[n-1]+phi_nco[n-1])/2 ) that's an oscillator at the mean frequency (omega[n] +G_nco[n]*x_nco[n])/2 which we hope to filter out with some LPF with transfer function H_lpf(z) = P + I*z/(z-1) + D*(z-1)/z which happens to also be our PID controller. P, I, and D are constants you might wanna fiddle with to get the response of this PLL servo-mechanism. P for "proportional", I for "integrator", and D for "differentiator". the second term sin((phi[n]-phi_nco[n])/2 is approx (phi[n]-phi_nco[n])/2 as the phase of the NCO closes in on the phase of the input, x[n]. recalling, phi_nco[n] = phi_nco[n-1] + G_nco*x_nco[n] that is an LTI system with Z transform Phi_nco(z) = (1/z)*Phi_nco(z) + G_nco*X_nco(z) so the transfer function of the NCO (input to phase) is H_nco(z) = Phi_nco(z)/X_nco(z) = G_nco/(1 - 1/z) = G_nco* z/(z-1) which is an integrator. so the input to our "plant" is (phi[n] - phi_nco[n])/2 or (Phi(z) - Phi_nco(z))/2 this is passed through the controller to get the input to the NCO X_nco(z) = H_lpf(z) * (Phi(z) - Phi_nco(z))/2 but Phi_nco(z) = H_nco(z) * X_nco(z). now, can you plug in the appropriate expressions and solve for the control system that processes phi[n] and results in phi_nco[n]? i'm starting to crap out. here, the "plant" already has an integrator in it. i think you'll find out that, if you crank it out, that the P constant will get attached to the in-loop integrator and the I becomes a 2nd-order integrator, and the D becomes the proportional term. and note that the phase discriminator is only that, you could dream up some algorithm (maybe using zero crossings) to get another term that is proportional to the frequency difference that would help your PLL close in on frequency differences when the frequencies are not very close. the principle value of the phase difference would be wildly varying and might get filtered out. but as the frequency difference is closed, then the phase difference is the only remaining non-zero term in the feedback. then this analysis kicks in. that's the best guess i have for how to think about a simple software PLL/ r b-j
fixing a couple mistakes that some proofing would have prevented...

__________________________________________________________

the NCO output:

    y_nco[n] = cos( phi_nco[n] )

where

    phi_nco[n] = phi_nco[n-1] +  G_nco*x_nco[n]

and x_nco[n] is the control "voltage" input for the NCO.  radians are
actually dimensionless, but G_nco is the "NCO gain" and might be in
units of radian/nco_ref or G_nco = 1/nco_ref.  when x_nco[n] =
nco_ref, the phase of the NCO advances 1 radian per sample.  but if
you consider x_nco[n] dimensionless, so is G_nco.

the input is

    x[n] = A * sin( phi[n] )

where

    phi[n] = phi[n-1] + omega[n]

for this simple analysis, omega[n] gotta be close to G_nco*x_nco[n].

so the output of the simple 4-quadrant multiplier is

    x[n]*y_nco[n] = A*sin(phi[n]) * cos(phi_nco[n])

    = A * sin(phi[n]) * cos(phi_nco[n])

    = A/2*( sin(phi[n]+phi_nco[n]) + sin(phi[n]-phi_nco[n]) )

the first term is

    A/2*sin(phi[n]+phi_nco[n]) = A/2 * sin(
                                 omega[n] + G_nco[n]*x_nco[n]
                             +   phi[n-1]+phi_nco[n-1]
                                 )

that's an oscillator at the sum of frequencies

    omega[n] + G_nco[n]*x_nco[n]

(about twice omega[n]) which we hope to filter out with some LPF with
transfer function

    H_lpf(z) = P + I*z/(z-1) + D*(z-1)/z

which happens to also be our PID controller.  P, I, and D are
constants you might wanna fiddle with to get the response of this PLL
servo-mechanism. P for "proportional", I for "integrator", and D for
"differentiator".

the second term

    A/2 * sin( phi[n] - phi_nco[n] )

survives the LPF and is approximately

    A/2 * (phi[n] - phi_nco[n])

as the phase of the NCO closes in on the phase of the input, x[n].

recalling,

    phi_nco[n] = phi_nco[n-1] +  G_nco*x_nco[n]

that is an LTI system with Z transform

    Phi_nco(z) = (1/z)*Phi_nco(z) + G_nco*X_nco(z)

so the transfer function of the NCO (input to phase) is

    H_nco(z)  = Phi_nco(z)/X_nco(z)

              = G_nco/(1 - 1/z) = G_nco* z/(z-1)

which is an integrator.

so the input to our "plant" is

    phi[n] - phi_nco[n]

which gets scaled by the phase discriminator to

    A/2 * (phi[n] - phi_nco[n])

or

    A/2 * (Phi(z) - Phi_nco(z))

this is passed through the controller to get the input to the NCO

    X_nco(z) = H_lpf(z) * A/2 * (Phi(z) - Phi_nco(z))

but

    Phi_nco(z) = H_nco(z) * X_nco(z).

now, you can plug in the appropriate expressions and solve for the
control system that processes phi[n] and results in phi_nco[n].

__________________________________________________________


just to add another note, i think the phase discriminator would need
to AGC the amplitude of the input, x[n] so that there would be no
variation of A and that would be a known constant that's specified
before the analysis.


r b-j
> >fixing a couple mistakes that some proofing would have prevented... > >__________________________________________________________ > >the NCO output: > > y_nco[n] = cos( phi_nco[n] ) > >where > > phi_nco[n] = phi_nco[n-1] + G_nco*x_nco[n] > >and x_nco[n] is the control "voltage" input for the NCO. radians are >actually dimensionless, but G_nco is the "NCO gain" and might be in >units of radian/nco_ref or G_nco = 1/nco_ref. when x_nco[n] = >nco_ref, the phase of the NCO advances 1 radian per sample. but if >you consider x_nco[n] dimensionless, so is G_nco. > >the input is > > x[n] = A * sin( phi[n] ) > >where > > phi[n] = phi[n-1] + omega[n] > >for this simple analysis, omega[n] gotta be close to G_nco*x_nco[n]. > >so the output of the simple 4-quadrant multiplier is > > x[n]*y_nco[n] = A*sin(phi[n]) * cos(phi_nco[n]) > > = A * sin(phi[n]) * cos(phi_nco[n]) > > = A/2*( sin(phi[n]+phi_nco[n]) + sin(phi[n]-phi_nco[n]) ) > >the first term is > > A/2*sin(phi[n]+phi_nco[n]) = A/2 * sin( > omega[n] + G_nco[n]*x_nco[n] > + phi[n-1]+phi_nco[n-1] > ) > >that's an oscillator at the sum of frequencies > > omega[n] + G_nco[n]*x_nco[n] > >(about twice omega[n]) which we hope to filter out with some LPF with >transfer function > > H_lpf(z) = P + I*z/(z-1) + D*(z-1)/z > >which happens to also be our PID controller. P, I, and D are >constants you might wanna fiddle with to get the response of this PLL >servo-mechanism. P for "proportional", I for "integrator", and D for >"differentiator". > >the second term > > A/2 * sin( phi[n] - phi_nco[n] ) > >survives the LPF and is approximately > > A/2 * (phi[n] - phi_nco[n]) > >as the phase of the NCO closes in on the phase of the input, x[n]. > >recalling, > > phi_nco[n] = phi_nco[n-1] + G_nco*x_nco[n] > >that is an LTI system with Z transform > > Phi_nco(z) = (1/z)*Phi_nco(z) + G_nco*X_nco(z) > >so the transfer function of the NCO (input to phase) is > > H_nco(z) = Phi_nco(z)/X_nco(z) > > = G_nco/(1 - 1/z) = G_nco* z/(z-1) > >which is an integrator. > >so the input to our "plant" is > > phi[n] - phi_nco[n] > >which gets scaled by the phase discriminator to > > A/2 * (phi[n] - phi_nco[n]) > >or > > A/2 * (Phi(z) - Phi_nco(z)) > >this is passed through the controller to get the input to the NCO > > X_nco(z) = H_lpf(z) * A/2 * (Phi(z) - Phi_nco(z)) > >but > > Phi_nco(z) = H_nco(z) * X_nco(z). > >now, you can plug in the appropriate expressions and solve for the >control system that processes phi[n] and results in phi_nco[n]. > >__________________________________________________________ > > >just to add another note, i think the phase discriminator would need >to AGC the amplitude of the input, x[n] so that there would be no >variation of A and that would be a known constant that's specified >before the analysis. > > >r b-j >
Thank you very much for this great explanation.