Math question - phase locking, filtering, etc.

Names have been changed to protect the innocent.

I have a spinning wheel in a dark room. On the circumference of the wheel are two lights spaced 180 degrees apart. At regular intervals, say 1/sec, one of the lights flashes (which one of the two is based on some signal). The flash is fast enough that I can’t track movement, but I can record the instantaneous position with a certain amount of noise (let’s say Gaussian). The wheel speed can change over time, and generally slowly, though there can be transients (like an initial wheel spinup).

The question is: how do I measure both the speed of the wheel and keep track of which of the lights flashes?

Some obvious points:

  • There’s no way to distinguish between light #1 and #2, since there’s no initial reference point, but other than the initial ambiguity I should be able to keep track.
  • The wheel can’t spin faster than 90 degrees/sec, else there would be no way to distinguish a flash moving in one direction vs. the opposite light moving in the other direction. So, no need to consider that case (in practice it’ll be much less).

My solution so far:
Keep a running estimate of the wheel position and angular velocity. When a new sample comes in, first update the predicted position with the predicted velocity. Take the difference in angle between the measured and predicted position. If it’s over 90 degrees in either direction, add/subtract 180 degrees and record that as a #2 flash; else don’t change the angle and record a #1 flash.

Correct the predicted position by averaging in the measured (corrected) position: say, Pp’ = Pp0.95 + Pm0.05. Compute a new predicted angular velocity similarly.
This works, though I suspect I can do better. I’m mainly fishing for ideas at the moment. One problem is that it’s hard to characterize the noise, and so picking the averaging constants is basically a random guess. I’d like something that is somehow more auto-adjusting. Ideas? Thanks!

If you had really sensitive instruments, you could calculate the redshift of the light, and figure out if it’s moving toward you or away from you. That would tell you if each flash is at the top or bottom, depending on how you’re situated relative to the wheel. I’m imagining you looking at the wheel edge-on.

As you realise, you have a phase locked loop. These are rather well understood beasts. Your only quirk is the use of 2 lights, but since you can disambiguate them, their presence doesn’t affect the working of the loop proper.

Each time you see a flash you work out where it was. Easy. Next estimate the error in the measured position relative to the predicted position. This is an error angle.

Now work out the error in the angular velocity. Now apply this error to the estimate of angular velocity (aka frequency) - but don’t do it all at once. You need to damp this update. A simple mechanism is to have an error term, and keep this updated, and every cycle you apply a portion of the error term to the frequency estimate. Subtract the amount of correction applied from the error term. This creates a very simple low pass filter that averages the error term over time. You can be more sophisticated than this. The actual rate at which the frequency estimate is updated is what is important.

What you need to do is choose damping parameters for the loop operation. The rate at which you change the estimate of frequency - you don’t want this to change instantly in each new estimate. That will couple the noise in position measurement straight into the position estimation. Choose a low pass filter time constant that filters out the position measurement noise as well as possible. Since you have some idea of the maximum rate of change of angular velocity you can estimate the time constant needed.

Since it is all in software, you can be more tricky.For instance you can start the system with a very wide filter, which is wide enough to cope with the startup transients. As the system comes up to speed and stabilises you can narrow down the filter.

To be a bit more precise in terminology.

When you compare the measured location with the predicted location - PLL theory calls this the phase detector.

Your estimator of angular position is a simple function of time and period. This is your oscillator. You can change the period, which makes it a variable oscillator.

The averaging of error update of the oscillator period is the loop filter.

I haven’t tried this - but this pagepreports to generate software PLLs directly from parameters. You could start here anyway.

Thanks for the tips!

Sounds pretty close to what I have now (just working things out from scratch). I’m semi-familiar with FIR and IIR filtering, and this kind of “y’ += a*(x - y)” filter is basically a trivial IIR filter. I can probably do a bit better… suggestions?

Definitely. An early version of my program did just this: the new position estimate was just the last measurement with the 180 degree correction. It worked, but very poorly. It did have infinite bandwidth, though (i.e., it would resync to a new wheel speed within one sample)!

For my problem, a glitch where one sample is wrong is not a big deal compared to a full, permanent phase flip. So coupling the noise into position is actually really bad.

Yeah, I think this is the approach I want to ultimately take. The problem is that the transients are fairly unpredictable, and so I’d like to infer them from the data. Maybe I can widen the filter if the positional error starts exploding?

Thanks much for the help here. I’ve read a bit about PLLs but am sometimes lost on the terminology and exactly how it maps to the problem at hand, so that was very useful.

Cool–I’ll check it out.

Sadly, there’s no actual wheel, or even lights :). My description is a fairly precise analogy for what I’m doing, but only in an idealized sense. The only inputs are a series of (x, y) measurements–there’s no additional information. In fact, I don’t even know the “diameter of the wheel”, and it might even change over time, but since angle is the only important variable here and its easily derivable from the coordinates, that’s not a problem.

10 virtual cookies for anyone that figures out what I’m actually doing :).

It has been a long time since I was in college and solved some of the problems like this.

Have you tried using the Simulink tool in Matlab ? It is a very intuitive tool to solve problems like this - you can even simulate the problem you are describing.

As to the solution - I would recommend a Neural Network that can learn and then predict the noise level in your measurement. It is fairly simple to implement such an ANN in Matlab. See this

I’ll take a look; thanks (if I can get my hands on Matlab). I’ll say that the solution has to be pretty fast–the real sample rate is on the order of a megahertz, and this part of the system can’t take more than several percent of a midrange CPU. So there are limits as to how fancy I can go.

I am not very abreast of the latest speeds in Data Acquisition Boards / Hard drives - in my time the fastest was a PCI board with a SCSI harddrive. I did go to pretty high sample rates though. Take a look at this DAQ board - it can do 10 MS/s.

If you can access Matlab it has toolboxes that are specifically for the purpose of designing control systems, and a PLL is a favourite. Sadly Matlab plus Simulink is not a cheap thing.

I use Scilab, which is an open source alternative to Matlab. It gets me by with most things I want to do. Google “scilab pll design” and you will turn up some useful stuff. Scilab runs on most platforms, and will run most Matlab code.

Mostly PLL design is directed at the analog domain, so op-amps and the like. But the actual parameters and theory behind doesn’t change. What you want to do is model the system to get you the right parameters.

This sounds like a hard real time system problem. That gets you into a whole new ballgame. A sample rate of MHz, is going to really stretch things unless you are very careful. You will be in the regime of needing guaranteed interrupt times, and so forth. Unless you have a dedicated RTS platform I don’t think is is going to be all that viable.

[spoiler]I’m working with a Software Defined Radio board which can acquire a pair of 12-bit samples (in-phase and quadrature) at 40 MS/s. Not bad for $400!

With my project, the signal in question will be frequency-shifted to DC, and the stream then decimated to ~1 MS/s, since I don’t need the full bandwidth. So the phase detector isn’t going at the full rate, although earlier stages in the pipeline may well be.[/spoiler]

Naw; there’s enough buffering that timing isn’t really an issue. I just need a reasonably low flop count per sample (tens, not thousands) to keep the throughput above a certain level. My current approach is just fine CPU-wise, and sounds pretty close to a baseline PLL, except without any theory backing the smoothing constants (i.e., I fiddled with them until they kinda worked). I can afford a few more cycles if a fancier filter proves necessary.

Thanks for the Scilab tip. I’ve used Octave a bit but the learning curve was pretty steep and the graphical side didn’t quite seem up to scratch.

I played around with PLLs and such and have some… results.

The main difference between my existing strategy and a “real” PLL, as far as I can tell, is in the frequency (angular velocity) estimation. I had been using an estimate of (currentSample - previousSample). That works, but has a few issues: first, that the noise is greater since I’m using two measured samples, and so the difference has more noise than either alone; and second, it averages in the phase shift (which looks like a really fast spike, and so throws off the average).

So I used the more standard approach of using (correctedSample - estimatedPhase) as the velocity, with correctedSample being currentSample with the 180 degree phase shift subtracted out if need be.

The problem is that the phase shift prediction only really works if it’s locked on already, and when it’s not locked on I get random prediction, which throws the PLL into chaos. So it just never reliably locks.

So I took a different approach. I knew the problem came down to having a better frequency estimator. The simple IRR filter worked ok, but not great–mainly because of those phase shifts acting like outliers and throwing things off. So what’s an outlier-resistant estimator? The median!

I found a technique for computing a moving-window median in O(log(n)) time per sample and used that instead, with a window of about a hundred samples. It works great! It locks on almost instantly and the drift is nearly zero.

Aside from the outlier rejection, this technique has the nice advantage that every sample within the window has equal influence, and every sample outside the window has zero influence. This is what I want since the presence of the signal is binary–either I’m getting good samples or not; there’s no notion of more recent samples being “better” than old ones.

Years later, I have run across a better answer.

It still requires a PLL, but with a preprocessing pass: treating the coordinates as a complex number, you square them.

This has the very convenient effect of making the two phases coincident. You can then run the PLL on that without the phase shifts throwing things off (this is very important if your bit rate is not much lower than the sample rate). It also has the effect of doubling the frequency, but that’s easy to get rid of later.

There are lots of ways of showing how it works, the most elegant probably being the Euler formula:
(e[sup]it[/sup])[sup]2[/sup] = e[sup]2it[/sup]
(e[sup]i(t + π)[/sup])[sup]2[/sup] = e[sup]2it + 2iπ[/sup] = e[sup]2it[/sup]e[sup]2iπ[/sup] = e[sup]2it[/sup]

It’s also trivial to work out with doing the complex arithmetic by hand, or via trig identities, but it all works out the same.

Another convenient effect: if I need to actually estimate a frequency, one way is to do a Fourier transform and look for the spike. But if the spike is at some fractional frequency, the accuracy isn’t that great. So instead I can do a few squarings to multiply the frequency by a power of two, then look for the spike. Dividing back out gives me a better fractional estimate.

The magic search term is carrier recovery.

Cool. And what real world project is/was this SDR PLL for? Or can’t you say?

Oh, sure–it’s actually pretty interesting (to me, at least).

Then, as now, it’s for a cubesat project. However, back in 2014 I was using it mostly for ground testing purposes and didn’t have stringent performance requirements.

Today, I’m working on an onboard receiver for an entry in the NASA Cube Quest Challenge.

The craft will reach fairly insane distances–tens of millions of kilometers out. And doesn’t even have a proper directional antenna. We will have time on the Deep Space Network or similar dishes for the Earth side, but we still have to pick a signal out from the noise on the craft side.

We’re using CCSDS-encoded packets, modulated in BPSK, and layered in a powerful error-correcting code.

BPSK stands for “binary phase shift keying”, and as you might imagine from previous posts, is a way of transmitting bits by transmitting one of two phases of a carrier wave, 180 degrees apart. It’s simple, reliable, and noise resistant.

However, as mentioned it’s a tad tricky on the receiving end. One of the first things you have to do with any demodulator is to synchronize with the transmitter (since the timings are never perfect–even quartz oscillators are roughly a part per million, and that’s 2 kHz when you’re at 2 GHz). But with BPSK, anything you do to synchronize with one phase tends to cancel out with the other phase, since they’re exactly opposite.

One approach is to detect when the phase shifts and alter how you update the PLL. But that starts to break down as the channel gets noisier and noisier.

What I just described is another approach–it transforms the waveform into one where both phases look the same. I run the PLL on that and I never have to discriminate between noise and bit transitions.

One minor flaw I’ve found: doubling the frequency can alias the higher frequencies (near the Nyquist limit) back into lower ones. If noise levels get really high, this can make the carrier harder to discriminate. The solution seems to be to apply a simple low-pass filter first to reduce that high-frequency noise energy.

All in all I’m pleased at how it’s coming along. I can inject noise at several times the signal level and still decode properly. By eye, you can’t even see that anything at all was transmitted, let alone see individual bits. Error correcting codes are pretty amazing things.

You have *much *better hobbies than I do. I’m envious, but far too lazy to act on it. Congrats on having so much more oomph. Thanks for sharing.

Ref the snip, there’s also the ever-changing Doppler shift between Earth & the cubesat. It’s a slow drift, but you can’t just tune up the nominal freq of the transmitter & expect to hear anything right there. Presumably you’ve got pretty solid position & velocity data so you can pre-compute the Doppler offset, but there’s noise in that calc too.

I notice your ref to CCSDS is pointing to an obsolete version of the spec. At the risk of being an amateur telling a pro his business, I sure hope that’s a deliberate version stabilization decision by somebody who controls the meat of the data chain.

Heh… I really don’t think of myself as having that much oomph, but somehow I have enough friends with that go-out-and-get-em attitude and they manage to rope me into their shenanigans. I think they know how to set nerd traps that will capture the problem-solving part of my brain. At any rate, it is pretty fun to work on this stuff.

Doppler is definitely a thing, but in principle the Deep Space Network has correction for that. Though this is an area where we need more information; their equipment seems largely designed around coherent operation–that is to say, the receiver can lock onto a carrier sent from Earth and retransmit it back. The Earth radio can then do very fine measurements of the frequency shift. I’m not sure we can get away with that mode of operation with our radio, so we may have to make do with less accurate ranging.

I could in principle just sample at a wide bandwidth and lock onto anything in a broad range. Doppler could add up to perhaps 45 kHz, so I just have to sample at somewhere above (twice) that rate to see the signal. Doable, but we also don’t have unlimited resources on the craft, so I’d rather avoid that.

It also affects the bit decoding–right now, I can decode a reasonably long packet just by stepping forward a bit at a time. The oscillator is good enough that it can stay in sync for a while. But with Doppler, the bits themselves are spaced differently, and I’d have to account for that.

It’s solvable but it’s still more complication. So it would be nice if the DSN correction did all that for us.

Good eye regarding the CCSDS doc; that just happened to be the first link I found, and didn’t notice the big “Historical Document” on the first page. We are in fact using a more recent version, though at quick glance not much has changed that affects us.

However, the DSN standards (which implement a subset of CCSDS) did change, and that ended up being pretty significant to us. The codes we thought we’d use were officially retired after the Cassini mission (i.e., under a year ago). They probably would have kept them working for a bit longer, but w’re flying on the SLS, and that’s been delayed, so now there’s no chance those codes would still be around. So we had to switch.

Lucky dog! We don’t have any available memory for a 100-point moving-window. And the computed PLL was working perfectly as long as the frequency was bounded… But there is a supplier out there who decided to use the frequency for signalling, and will run the frequency at the limits provided by the standard for error conditions! I think we are adding extra board components to do zero-point-crossing detection now.

I know, right? Although I do enjoy to-the-metal optimization, shaving off every cycle and byte of memory, sometimes it’s nice to be able to just implement the “right thing” and not have to worry about HW constraints.

With the latest project, I’m doing a 262,144-sample floating point FFT in to pick up the carrier. No biggie; I’ve got a few hundred MB of RAM and a quad-core ARM processor with an FPU.

Zero-crossing detectors are great if the noise isn’t too bad. If the data has weak-to-no error correction, that kind of detector is almost certainly good enough, since if it were worse you’d be getting junk data out of it anyway. Not an option for the current project, though.