Coherent sampling is quite difficult to meet under the lab conditions. One has to go for windowing to characterize the dynamic performance of ADC. Though seminal papers and reports [1-2] lie on my desk for quite long time, I still feel difficult to really understand windowing. Recently, I found two posts from Nerd Rage (Windows of Opportunity & ENBW) quite helpful, which are easier to digest.

And now, finally, I decide to write down my very limited understanding on FFT with windowing, focusing on characterizing the dynamic performance of ADC in Matlab. Let’s go straightforward by first looking at the Matlab script I am using to calculate the SNR/SNDR of ADC.

fs = 1e6;             % Sample rate
nfft = 1024;          % Number of FFT
cycles = 113;         % Number of input periods
fin = cycles/nfft*fs; % Input frequency
data = data(1:nfft);  % Number of data to be the same as nfft

w = 0.5*(1 - cos(2*pi*(0:nfft-1)/nfft)); % Hann window
cg = sum(w)/nfft;                        % Normalized coherent gain
enbw = sum(w.*w)/(sum(w)^2)*nfft;        % Normalized equivalent noise bandwidth
nb = 3;                                  % Signal bins
dcbin = (nb+1)/2;                        % Number of DC bins

if (size(data,1)~= size(w,1))            % Check dimention
   w = w';
end

ss = abs(fft(data.*w));      % FFT with windowing
ss = ss/nfft/cg;             % Compensate for window attenuation
ss = ss(1:nfft/2).*2;        % Drop the redundant half but keep total power the same

signal_bin = cycles+1;                           % Signal bin, Matlab array starts from 1
dc_bins = 1:dcbin;                               % DC bins
all_bins = setdiff(1:nfft/2, dc_bins);           % Disregard DC bins
signal_bins = signal_bin + (-(nb-1)/2:(nb-1)/2); % Signal leakage bins
other_bins = setdiff(all_bins, signal_bins);     % Further discard signal bins

fh = (2:10)*fin/fs;                              % Harmonic tone: (-/+m)fin + (-/+k)fs
   while max(fh) > 1/2
   fh = abs(fh - (fh > 1/2));                    % If harmonic tone fh>fs/2, it is equal to fs-fh
end
harm_bins = round(nfft * fh) + 1;                % Harmonic bins (2nd - 10th)
harm_binsl = zeros(length(harm_bins),nb);        % Find Harmonic leakage bins
for i = 1:length(harm_bins)
    harm_binsl(i,:) = ((harm_bins(i) + (-(nb-1)/2 : (nb-1)/2)));
end
harm_binsl=reshape(harm_binsl',length(harm_bins)*nb,1);  % Convert matrix to array
harm_binsl=unique(harm_binsl);                           % Discard the repetitive harmonic bins

noise_bins = setdiff(other_bins, harm_binsl);            % Further discard the harmonic bins

Psignal = sum(ss(signal_bin).^2);                        % Signal power
PnoiseD = sum(ss(noise_bins).^2)/enbw/length(noise_bins);% Noise PSD
Pnoise = PnoiseD*length(all_bins);                       % Total noise power
Pharm = sum(ss(harm_bins).^2);                           % Power of harmonics

snr = 10*log10(Psignal/Pnoise);                          % Calculate SNR
sndr = 10*log10(Psignal/(Pnoise+Pharm));                 % Calculate SNDR
enob = (sndr - 1.76)/6.02;                               % Calculate ENOB
thd = -10*log10(Psignal/Pharm);                          % Calculate THD

Pharm_max = max(ss(harm_bins)).^2;
Pnoise_max = max(ss(noise_bins)).^2;
sfdr = 10*log10(Psignal/max(Pharm_max,Pnoise_max));      % Calculate SFDR

sdb = 10*log10(ss.^2);
f = (1:length(ss))/nfft;                                 % Frequency vector normalized to fs

plot(f, sdb, 'k-','linewidth',1.5);
xlabel('Frequency [ f / f_s ]','FontSize',10);
ylabel('Power Spectrum [ dB ]','FontSize',10);
grid on;
text(0.3, -40,...
sprintf('SNDR = %.1fdB \n SNR = %.1fdB \n THD = %.1fdB \n SFDR = %.1fdB \n ENOB = %.1fbits', sndr, snr, thd, sfdr, enob), ...
'FontSize',10);

An example of FFT plot will look like this:

Fig.1 An example of FFT plot showing the dynamic performance of ADC

A default FFT using rectangular/box window has a normalized coherent gain as 1, a normalized equivalent noise bandwidth as 1, and only 1 signal bin. Hence, if other types of window are used, the corresponding window properties need to be specified. It’s all about two values: 1) the sum of the window terms; and 2) the sum of the squares of the window terms. Equations from the seminal paper [1] tell us why the two sums matter.

Let the input sampled sequence be defined by

where q(nT) is a white-noise sequence with variance . Then the signal component of  the windowed spectrum is given by

Hence, the output amplitude of the noiseless signal is the input amplitude multiplied by a term which is the sum of the window terms (S1). This term is called the processing gain (sometimes called coherent gain) of the window. The rectangular window has the largest gain compared to other windows. In the Matlab script, the normalized coherent gain (normalized by S1 of the rectangular window) is specified.

The noise component of the windowed spectrum is given by

The noise power is calculated using the expectation operator

As additive noise is assumed to be white, the above value represents the noise floor level (or the noise power spectral density), which is also constant. Notice the power gain of noise is the sum of the squares of the window terms (S2).

The noise bandwidth is calculated by the total noise power (sigma^2*S2*fs) divided by the peak power gain of the window (S1^2). Here we only focus on the multiplication term to the input noise (S2*fs/S1^2). Further introducing a parameter, fres, the width of one frequency bin (fs/N), then the normalized equivalent noise bandwidth (ENBW) is given by

Therefore, to obtain correct power levels of the signal and the noise, the normalized coherent gain is subtracted from the PSD and the normalized ENBW is subtracted from the calculated total noise power, respectively.

Nerd Rage explains the ENBW from window energy point of view, which is quite insightful. Allow me to copy one of his picture and some remarks here:

Fig. 2 Some properties of the Hann window and its instrument function. The window function (top left) is normalized to T = 1; graph also shows the window energy relative to the box window. Instrument function (top right) is shown for the spectral region covering the first 10 DFT bins. Close-up (bottom) shows main lobe height, maximum scalloping loss (see below), relative height of first sidelobe, and fraction of energy contained in main lobe, E0, to total energy E. (Courtesy: Nerd Rage)

“For the Hann window the main lobe height is -6.02dB and therefore the height of any single spectral line will be 6.02dB below its real value. To obtain correct energy levels of the spectral peaks (in the absence of scalloping), the main lobe height (in dB units) is usually subtracted from the PSD. However, this overcompensates the window-specific reduction of the noise floor – for the Hann window, peak compensation is 6.02dB while the noise floor is only 4.26dB below its true value. This decreases the observed SNR by 1.76dB or a factor of 1.5.”

Note: for Hann window, 20log(s1/N) = -6.02dB; 10log(s2/N)=-4.26dB; N*s2/s1^2 = 1.5.

Regarding the number of spread bins for a certain window function, a simple Matlab script performing FFT with windowing can do the job:

% Here you create an input 
fs = 1024;
nfft = 1024;
fin = 16;
cycles = 16;
vin = sin(2.*pi.*fin.*(0:1:nfft-1)./fs);
% Here you perform FFT with and without windowing
ss = abs(fft(vin(1:nfft)))/(nfft/2);
w = 0.5*(1 - cos(2*pi*(0:nfft-1)/nfft));
cg = sum(w)/nfft;
sshann = abs(fft(vin(1:nfft).*w))/(nfft/2)/cg;
% Here you do the plot
figure(1)
subplot(1,2,1)
stem(ss);
title('Rectangular');
xlim([1 25]);
grid on;
subplot(1,2,2)
stem(sshann);
xlim([1 25]);
title('Hann');
grid on;

Then it comes the plot:

Fig.3 FFT with Rectangular and Hann window (For Hann window, signal bins = 3)

Sorry for the boring and not well-written post (though I’ve tried my best). Fig.4 is just for a smile when one finishes reading all the above dull words.

Last but not least,  may Mr. Fourier’s wisdom be always with us!

Fig.4 A sketch of Mr. Fourier

[1] Frederick J. Harris, “On the Use of Windows for Harmonic Analysis with the Discrete Fourier Transform,” Proceedings of the IEEE, vol 66, pp. 51-83, January 1978.

[2] G. Heinzel, A.Rudiger, and R. Schilling, “Spectrum and spectral density estimation by the Discrete Fourier transform (DFT), including a comprehensive list of window functions and some new flat-top windows,” February 15, 2002.