Doppler Spread: It's the Small Things

Overview

Rob Robinett, AI6VN, recently reached out to me to let me know that he was seeing Doppler spreads of zero for some signals using the wsprd.c patch I wrote. I’d never seen that before, so Rob sent me a WAV file with some very, very clean signals that yielded a Doppler spread of zero. After taking a look at how they were processed, I learned a little bit more about why you have to pay particular attention to the FST4 algorithm for small Doppler spreads.

Overshoot Correction

Let’s take a look at a made-up, simplified FFT of g (channel gain, signal times the complex conjugate of a reference signal) for a signal. Remember that the FFT (whether it comes from FFTW or Numpy) is an array where each index represents a bin of frequencies:

fft

To calculate this spread we’d start by taking the total power, 15 (1 + 2 + 9 + 2 + 1), and find the bins where our sum reaches 25% of 15, 3.75, and 75% of 15, 11.25.

Here’s a Python snippet:

total_power = 15
i25 = None
i75 = None
curr_power = 0
for i in range(0, 5):
    curr_power += fft[i]
    if curr_power >= 0.25 * total_power and i25 == None:
	i25 = i
    if curr_power >= 0.75 * total_power and i75 == None:
        i75 = i

If we ran this on our example we’d get i25 of 2 and an i75 of 2. Subtract them and you’d get zero. Multiply that difference by the frequency space each bin represents, DF (more on this number in Small Value Fix), and you’d still get zero!

So what’s going on here? Well our algorithm is overshooting where the power actually reaches its mark because it’s quantized by the integer of the bin. What we really want is an estimation of where the power reaches its mark in between the bins. FST4 does this through linear interpolation of the data:

interpolation

So if we assume the power rises linearly between the bins we can generating a floating point number that does a better job of representing where we reach our power mark. FST4 does this with the following equation:

$$OvershootCorrection = {(ExpectedPower - PreviousPower) \over (CurrentPower - PreviousPower)}$$

It should be noted that as of the writing of this article, FST4 was using a slightly incorrect equation for this. It will be fixed in the next version.

Now we can account for the overshoot correction from the previous index, remember our index travels past where we actually want to stop, and get a floating point number that should prevent our difference from being zero. Let’s add that to our example Python code and see what happens:

total_power = 15
x25 = None (1)
x75 = None
curr_power = 0
prev_power = 0 (2)
for i in range(0, 5):
    curr_power += fft[i]
    if curr_power >= 0.25 * total_power and i25 == None:
	x25 = i - 1 + (curr_power - 0.25 * prev_power) / (curr_power - prev_power) (3)
    if curr_power >= 0.75 * total_power and i75 == None:
        x75 = i - 1 + (curr_power - 0.75 * prev_power) / (curr_power - prev_power)
    prev_power = curr_power
1 Notice I changed the variable prefix to x as i being used for a non-integer pains me!
2 We also need to track the previous power. This may not work as expected if g reaches 25% power in the first bin. In that case, it might be reasonable to assume the signal spread out of the power region.
3 Here is the overshoot correction

If we ran this for 25%, we would get an overshoot correction of:

Overshoot Correction = (3.25 - 3) / (12 - 3) = 0.08
x25 = 1.08

And for 75%, we would get an overshoot correction of:

Overshoot Correction = (11.25 - 3) / (12 - 3) = 0.92
x75 = 1.92

This yields a difference of 0.84 which multiplied by DF will give us a non-zero Doppler Spread!

Small Value Fix

FST4 has one more trick up it’s sleeve regarding small Doppler Spread. It uses the following equation (in two parts):

$$w50 = \sqrt{1+diff^2}\times{df}$$

The comment next to it states, "Keep small values from fluctuating too widely."

Recall that w50 is another name for Doppler Spread and df is the range of frequencies contained within each FFT bin. df can be calculated as the sample rate divided by the amount of samples, which for wsprd.c is 375/46080 and for FST4 is 12000/1440000. This equation will clamp down Doppler Spread values as they approach zero, such that the smallest Doppler Spread value will be df. So for WSPR, don’t expect a Doppler Spread smaller than 0.008.

Conclusion

There’s one footnote in my FST4 Doppler Spread Algorithm in WSJT-X article regarding overshoot correction and it reads:

I found this doesn’t make a huge difference and I left it out of my Python implementation.

Boy was I wrong! I’ve added overshoot correction and the small value fix to my WSPR patch and wspr_spread.py scripts.