Spectral analysis with high selectivity but low number of points

I'm trying to create spectral plots similar to one below (it's from the datasheet of an AD1986 sample rate converter, page 8): In the comments below the plot, it says "Wideband FFT plot, 16K points".

I'm particularly interested in having high spectral selectivity (if that's the right term), with the 1kHz peak being very narrow, with 'only' having 16K FFT points.

I'm doing my exercise on a 1000.621 Hz sine wave, quantized to 16 bits with a 44.1kHz sample rate.

This code to create all the plots below can be found here.

Step 1: straightforward FFT of 1s of samples, with and without window: There are 2 issues with this:

• there's a lot of variability in the noise
• there's still a lot of spreading around the peak, despite using a Blackman window

Step 2: Average 100 FFTs over 100s of samples This averages the noise, and it overlays perfectly with the original FFT in terms of signal power.

Step 3: Do a single FFT over 100s of samples Increasing the number of samples by 100, has reduced the width of the FFT bins.

Of course, the noise isn't averaged, so I'm back to square one on that one. What's more, the noise floor due to quantization will be lower than for step 2, due to the increased number of bins.

Step 4: Merge 4410000 FFT bins of step 3 back to 44100 FFT bins

I could average a bunch of these to reduce the noise, just like in step 2, but the AD1896 image says it's a 16K point FFT. Which means that they used a different method.

So I think that I should also be able to merge 100 FFT bins back to 1 bin, so that I get the combination of both a narrow peak and a noise floor that that is similar in value to the one in step 2.

Here's the result after adding the absolute values of each bucket of 100 FFTs bin into 1 FFT bins: This is starting to look good, but the peak value is now well above 0dB. And the noise floor around ~-128dB is much higher than the one that was around ~-145dB of step 2.

What I really want is the result of step 2 with the same level of noise floor, but without the wide peak of step 3.

I have the feeling that adding the absolute values of each bin is not the right approach, but I don't know which would be the right approach, and why?

The best way to do ae this measurement: Use a sine wave generator that is phase locked to your data acquisition clock with a frequency that's an integer multiple of your sample rate divided by the FFT length. If you do this, the period of your sine wave becomes an integer number of samples, you don't get any spectral spreading and you don't need any windowing at all.

You combine different frequency bins by summing the energy of bins.

• About the integer number of periods: that's a suggestion that I've seen being made in documents about ADC performance testing. See section 2.2 of this document: mit.edu/~klund/A2Dtesting.pdf, for example. One thing is explicitly mention, is that this integer multiple should be odd and prime, to avoid there being a common number between the number of input period and the number of sample points. Feb 22 '21 at 17:21
• I've also seen the exact opposite recommendation being made: explicitly avoid an integer multiple. So it depends on that kind of measurement you want to do. It's not 100% clear to me when to do the integer multiple and when not. In this particular case, since I want to eventually reproduce the measurement results of the AD1986, where source and sink have asynchronous clocks, I think it's impossible to satisfy the multiple integer rule on both sides. Feb 22 '21 at 17:25

While I'd still like to know the answer to my question about how to merge multiple FFT bins into one, the real solution to increasing the selectivity was much easier to solve: I just had to use a different windowing function.

Here's the result of step 2 for 3 cases: no windowing function (blue), using the Blackmann window that I used earlier (green), and using a Kaiser window with a beta factor of 20 (orange): In my numpy code, the change is as simple as this:

w = np.blackman(len(y)/num_ffts); corr = len(w)/sum(w)

by this:

w = np.kaiser(len(y)/num_ffts, 20); corr = len(w)/sum(w)
• FWIW, I’m with you on this one, and I think you’d be within your rights to contact ADI or post on their support forum to find out how they generated their plots. Their plot at a glance seems to have fewer than 16k points, and the width of the peak is so narrow, their either using some slick windowing or they matched the input frequency and FFT length to get it dang close to a bin. Feb 23 '21 at 3:33