[time-nuts] A simple sampling DMTD

[Joseph Gwinn]

JD,

On Thu, 28 Nov 2019 12:00:01 -0500, time-nuts-request at lists.febo.com 
wrote:
> Send time-nuts mailing list submissions to
> 	time-nuts at 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 at lists.febo.com
> 
> You can reach the person managing the list at
> 	time-nuts-owner at 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 at gmail.com>
> To: Discussion of precise time and frequency measurement
> 	<time-nuts at lists.febo.com>
> Subject: Re: [time-nuts] A simple sampling DMTD
> Message-ID:
> 	<CAEoGJMnToj26awb=1FKFtJuGSduDtGB_2XBwFYx6D72bJwAEZw at mail.gmail.com>
> Content-Type: text/plain; charset="UTF-8"
> 
> Dear Joe,
> 
> On Wed, Nov 27, 2019 at 7:22 PM Joseph Gwinn <joegwinn at 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
> ******************************************