Re: [time-nuts] A simple sampling DMTD
JD,
On Thu, 28 Nov 2019 12:00:01 -0500, time-nuts-request@lists.febo.com
wrote:
Send time-nuts mailing list submissions to
time-nuts@lists.febo.com
To subscribe or unsubscribe via the World Wide Web, visit
http://lists.febo.com/mailman/listinfo/time-nuts_lists.febo.com
or, via email, send a message with subject or body 'help' to
time-nuts-request@lists.febo.com
You can reach the person managing the list at
time-nuts-owner@lists.febo.com
When replying, please edit your Subject line so it is more specific
than "Re: Contents of time-nuts digest..."
Today's Topics:
1. Re: A simple sampling DMTD (Jan-Derk Bakker)
Message: 1
Date: Wed, 27 Nov 2019 23:46:32 +0100
From: Jan-Derk Bakker jdbakker@gmail.com
To: Discussion of precise time and frequency measurement
time-nuts@lists.febo.com
Subject: Re: [time-nuts] A simple sampling DMTD
Message-ID:
CAEoGJMnToj26awb=1FKFtJuGSduDtGB_2XBwFYx6D72bJwAEZw@mail.gmail.com
Content-Type: text/plain; charset="UTF-8"
Dear Joe,
On Wed, Nov 27, 2019 at 7:22 PM Joseph Gwinn joegwinn@comcast.net wrote:
[snip]
The 1001-point FIR band pass filter is a good reference to get an idea of
the best case performance of the system, but it is computationally
infeasible to run on an 8-bit processor. While a cheap comb filter can
take
a bite out of the HF noise, canceling offset/drift is harder. Early on I
was looking into ways to average all samples of a single period of the
beat
note, but I had trouble finding a closed form solution to the fact that
the
beat note period is never an exact integer multiple of the sampling rate.
Through numerical methods I did end up finding a good estimator for the
"leakthrough" caused by the fractional part of the beat note period (as a
function of the measured period and the phase offset), which was fairly
inexpensive to implement. [snip]
In the radar world, the standard solution to the leakthrough problem is
to batch the data and apply a windowing function to the data in each
batch. Typically, the batches overlap such that every sample appears
in two batches. The window functions largely eliminate the splice
error due to the FFT, which in fact splices each batch into a circle,
causing a discontinuity at the splice. If the window function is very
small at the splice, there is little discontinuity power sprayed into
innocent FFT bins.
This is an interesting approach, and one which hadn't occurred to me while
I was working on the offset/drift cancellation. After some consideration I
am afraid it will not help much in this particular case.
What I want is to filter out as much DC/LF energy before running the zero
crossing detector (ZCD) on the rising edge of the sampled sine wave. To do
this right now I also run the ZCD on falling edges, and I remember the
(sub)sample position of the current and previous falling edges. I then
average a full nominal period of the sine wave between these two falling
edges(200 samples at 10Hz beat frequency and 2ksps after the first sinc^2
filter/decimator), to determine the offset at the rising edge in between.
This is indeed mathematically equivalent to calculating bin 0 of a forward
Fourier transform with a rectangular window. Because the actual period is
not an exact integer number of samples, in practice this bin is
contaminated by spectral leakage. Windowing would help here, if it weren't
for the fact that all of my wanted signal energy is in bin 1. As any
windowing function other than the rectangular window trades side-lobe
suppression for main lobe width, the cure would likely be worse than the
disease. This could be ameliorated by averaging over multiple periods
(increasing the transform length), but that would necessarily spread the LF
noise energy over multiple bins, again making it all but impossible to get
optimal rejection.
What I've done instead is make an estimator for the leakage. For a period
of 200 +/- 1 sample, the following has an estimator error of better than
-40dB:
alpha = 0.5*perOffset^2 + (1.5 - phase)*perOffset - phase
correction = (1-alpha) * pointN + alpha * pointN1
where
correction is the estimator (in sample values)
perOffset is the difference between the measured period (199...201
samples) and the nominal period (200 samples)
phase is the estimated starting phase of the averaging interval, in
fractional samples (range: [0,1>)
pointN is the value of the sample immediately after the averaging
interval (ie sample[200], when the averaging interval is
sum(sample[0]...sample[199])
pointN1 is the value of the second sample after the averaging interval
(ie sample[201]).
This works very well, both on simulated and actual data. I am embarrassed
to admit that I cannot explain why it must be this particular formula
(and I'm very open to get hints here).
Sincerely,
JD 'if you can't do the math, run simulations and do curve fits until it
works' B.
I'll have to think about it, but I will mention that with batch
processing using window functions, it's common to precede the FFT using
a simple FIR filter to eliminate the low-frequency energy (due to
clutter, DC leakage/offset et al), the problem being that the FFT alone
may not have sufficient rejection to prevent low-frequency energy
breakthrough.
An alternate approach is to compute the mean of the windowed data and
subtract that mean from all samples in the window before computing the
FFT. In analog terms, this is like a big coupling capacitor blocking
DC.
Joe G
End of time-nuts Digest, Vol 184, Issue 37
Hi,
On 2019-11-28 21:18, Joseph Gwinn wrote:
JD,
I'll have to think about it, but I will mention that with batch
processing using window functions, it's common to precede the FFT using
a simple FIR filter to eliminate the low-frequency energy (due to
clutter, DC leakage/offset et al), the problem being that the FFT alone
may not have sufficient rejection to prevent low-frequency energy
breakthrough.
An alternate approach is to compute the mean of the windowed data and
subtract that mean from all samples in the window before computing the
FFT. In analog terms, this is like a big coupling capacitor blocking
DC.
You also wants to remove the slope between the two sample end-points, or
else that slope represents an in-phase sawtooth function. A
window-function tends to do that.
Cheers,
Magnus
Dear Magnus,
Removing the slope between the two sample end points (or: trimming/adding
the fractional sample part of the period) is the point of the estimator I
posted earlier (
http://lists.febo.com/pipermail/time-nuts_lists.febo.com/2019-November/098450.html
).
In general: as much as I like having it in my toolbox, I don't see how
using an FFT would be the best tool for the job in a zero-crossing detector
for a DMTD, let alone this particular sampling DMTD. For one, this 8-bit
processor doesn't have the spare cycles to run FFTs on the 32-bit data I
get from my CIC^2 decimator; besides that, I would only be interested in a
single bin (the beat frequency), where it would be more efficient to simply
I/Q-demodulate the samples in software (O(N) vs O(N log N)). While I admit
that in the latter case windowing would help, at this point I/Q
demodulating (effectively calculating only a single bin of the DFT) does
not appear to have advantages over least squares fitting the arcsine of the
incoming samples. Am I missing something here?
Sincerely,
JDB.
On Fri, Nov 29, 2019 at 2:32 AM Magnus Danielson magnus@rubidium.se wrote:
Hi,
On 2019-11-28 21:18, Joseph Gwinn wrote:
JD,
I'll have to think about it, but I will mention that with batch
processing using window functions, it's common to precede the FFT using
a simple FIR filter to eliminate the low-frequency energy (due to
clutter, DC leakage/offset et al), the problem being that the FFT alone
may not have sufficient rejection to prevent low-frequency energy
breakthrough.
An alternate approach is to compute the mean of the windowed data and
subtract that mean from all samples in the window before computing the
FFT. In analog terms, this is like a big coupling capacitor blocking
DC.
You also wants to remove the slope between the two sample end-points, or
else that slope represents an in-phase sawtooth function. A
window-function tends to do that.
Cheers,
Magnus
time-nuts mailing list -- time-nuts@lists.febo.com
To unsubscribe, go to
http://lists.febo.com/mailman/listinfo/time-nuts_lists.febo.com
and follow the instructions there.
Am 29.11.19 um 11:45 schrieb Jan-Derk Bakker:
In general: as much as I like having it in my toolbox, I don't see how
using an FFT would be the best tool for the job in a zero-crossing detector
for a DMTD, let alone this particular sampling DMTD. For one, this 8-bit
processor doesn't have the spare cycles to run FFTs on the 32-bit data I
get from my CIC^2 decimator; besides that, I would only be interested in a
single bin (the beat frequency), where it would be more efficient to simply
I/Q-demodulate the samples in software (O(N) vs O(N log N)). While I admit
that in the latter case windowing would help, at this point I/Q
demodulating (effectively calculating only a single bin of the DFT) does
not appear to have advantages over least squares fitting the arcsine of the
incoming samples. Am I missing something here?
I admit that I did not follow this thread closely, but the Goerzel filter
is the single output line DFT , with O(n).
<
https://www.mathworks.com/help/signal/ref/goertzel.html;jsessionid=8816f77eb76ad7dd1913c7021698
< https://en.wikipedia.org/wiki/Goertzel_algorithm >
If you need to simulate floats, fractional integers are easiest.
I/Q demodulation probably requires to recreate a clean carrier if you want
absolute phase and not only relative jumps. That sounds more like FPGA
than 8051, or whatever the 8 bit processor of the day may be.
regards, Gerhard
OK, in a previous life I did build a system for geophysics, where they fed
dangerous amounts of AC into the soil and measured the potential at
some 50 nodes. Rubber boots required.
Each node had a 8951 to control some switches and communicate setups.
They were most exited when I gave them the sources so they could
implement FFT pre-processing locally on each node themselves.
That required willingness to suffer.