[SOLVED] SNR of a Quantized Signal using FFT

[iffe]

Hi all,

I have a problem in the calculation of SNR.

Sine wave is quantized and then I have calculated its SNR using FFT. The coherence condition with windowing as well is implemented to avoid leakage in the adjacent frequency bins.

I am unable to calculate the SNR accurately with different resolution of quantizer (ADC).

If I use a standard formula for quantization noise power (Q^2 / 12), the SNR corresponds with the theoretical formula of SNR (SNR = 6.02*NoOfBits + 1.76) but when I try to use my own calculated noise power, the SNR does not corresponds with the theoretical formula with different resolution.

I think, I am not calculating the noise power in the correct way. I have considered the noise spectrum starting from 2nd bin to input bandwidth.

By the way I am using the formula of SNR = 10*log10(RMS_of_Signal/Mean_of_RSS_of_Noise)

Please help.

/Iffe

 

[JoannesPaulus]

The signal power should be calculated as:

sig_pwr=sum(X(bin_start:bin_end).*conj(X(bin_start:bin_end)));

where X is the DFT of your signal, bin_start is the first bin of your DFT containing signal power, bin_end is the last bin containing signal power (for example, if you use a window that spreads the signal power across 7 bins and the center frequency is at, say, bin 191, bin_start=188; bin_end=194.
The noise power is calculated the same way but skipping the signal bins. Moreover, if you have a non-linear system, you should also skip the bins containing harmonic power (unless you are interested in SNDR).

 

[iffe]

Hi JoannesPaulus,
Many thanks for your kind reply.
I think I am doing by the same way.

I am using "hann" window and it is giving me signal power in the three bins of FFT. I have scaled the FFT magnitude by dividing with N/2 (the value of the FFT is determined in this way) in order to compensate with the window effect. Also, it seems that power of signal is calculated in the right way and that is why when I use ideal quantization noise (Q^2/12), it gives me correct results of SNR. I have neglected the main lobe bins while calculating noise power.

My simulation code is attached here with. Kindly have a look and do let me know if you can get some time.

%*
close all;% closing all the open figures
clear all;% clearing all previous variables
% =
% Generate a simulation signal
% =
fin = 257;
fs = 8192;
N = 8192;
adc_resolution = 16; % Or Quantizer resolution

Nw = floor((fin*N)/fs);
amp = 0.5;
signal = amp*sin(2*pi*Nw/N*[0:N-1]);

figure(1);
subplot (2,1,1);
plot(signal);
title('Original Signal');
ylabel('Signal Amplitude');
xlabel('Time');
%=
%%
% =
% Scalar quantization is implemented over here
% =
no_quantiz_levels = power(2, adc_resolution);
%As signal contains maximum value of 1 and minimum value of -1
partition = linspace(-1, 1, no_quantiz_levels);
codebook = linspace(-1, 1, no_quantiz_levels);
index = quantiz(signal, partition);

quantize_signal = zeros(1, length(index));
for i = 1:1:length(index)
quantize_signal(i) = codebook(index(i) + 1);
end

subplot (2,1,2);
stem(signal);
hold on
plot(quantize_signal, 'g*');
title('Original Signal and Quantized Signal');
legend('Original','Quantized');
% =
%%
% =
% FFT and SNR measurement
% =
%Remember that hann window reduces the amplitude of the signal to half in
%frequency domain and that is why scaling is done by N/4 instead of N/2
fft_signal = fft(quantize_signal.*hann(N)',N)/(N/4);
freq_scal = linspace (0, 0.5, N/2); % for normalized frequency

figure(2);
plot(freq_scal, abs(fft_signal(1:N/2)), '-x');
hold on
plot(freq_scal(Nw:Nw+2), abs(fft_signal(Nw:Nw+2)), 'r*');
title('FFT Plot (N = 8192, Hann window is used)');
xlabel('Frequency (Normalized)');

signal_value = ((sum((abs(fft_signal(Nw:Nw+2))).^2))/3)

noise_bins = [fft_signal(2:fin-1) fft_signal(fin+3:N/2)];
noise = (sum((abs(noise_bins)).^2)) / ((N/2)-3)

noise_formula = ((1/power(2,adc_resolution))^2)/12

SNR = 10*log10(signal_value/noise)
SNR_byQuantizationFormula = 10*log10(signal_value/noise_formula)

You can easily see that if I change the resolution the "SNR_byQuantizationFormula" gives the exact value while other gives not.
I am anxiously waiting for your reply.

/Iffe

 

[JoannesPaulus]

I noticed several problems. The quantizer is not correct, a simpler way to implement it is using the floor function:

quantize_signal=floor((no_quantiz_levels-1)*signal)+no_quantiz_levels/2

This will generate a quantized signal between [0 255].
You can calculate the scaling factor due to windowing by:

window_type=hann(N)';
window_scaling_factor=sum(window_type);

with this, you can use any window.
Moreover,

signal_value = ((sum((abs(fft_signal(Nw:Nw+2))).^2))/3)

Why are you dividing the signal power by 3?

I do not have matlab with me, I will try to run your code tonight...

 

[JoannesPaulus]

I also noticed that you are only taking half of the signal power and half of the noise power (two errors that cancel each other so, no harm there).

 

[iffe]

Hi Joannes Paulus,

Thanks for your email and your precise attention on my code. I have tried to implement your quantizer but it is not working correctly with my code (some problem is with my understanding). If you run my code, I hope you will note that my quantizer is working well. I think that the window scaling factor and Signal power calculation is correct as "SNR_byQuantizationFormula" gives a correct value while using ideal quantization noise power.
I think it will be better if you can send me your recommended changes in a Matlab file so that I can see the difference between your and mine code.


With Regards:
Ahsan

I also noticed that you are only taking half of the signal power and half of the noise power (two errors that cancel each other so, no harm there).

 

[iffe]

Hi Joannes Paulus,

I have added the changes according to your requirement but still the results are the same and also the SNR calculated with ideal quantization power does not corresponds to theoretical formula.
The code is attached here with.

%*
close all;% closing all the open figures
clear all;% clearing all previous variables
%=
% Generate a simulation signal
%=
fin = 257;
fs = 8192;
N = 8192;
adc_resolution = 16; % Or Quantizer resolution

Nw = floor((fin*N)/fs);
amp = 0.5;
signal = amp*sin(2*pi*Nw/N*[0:N-1]);

figure(1);
subplot (2,1,1);
plot(signal);
title('Original Signal');
ylabel('Signal Amplitude');
xlabel('Time');
%=
%%
%=
% Scalar quantization is implemented over here
%=
no_quantiz_levels = power(2, adc_resolution);
quantize_signal = zeros(1, length(signal));

quantize_signal = floor((no_quantiz_levels-1)*signal);%+no_quantiz_levels/2

subplot (2,1,2);
stem(signal);
hold on
plot(quantize_signal, 'g*');
title('Original Signal and Quantized Signal');
legend('Original','Quantized');
%=
%=
% FFT and SNR measurement
%=
%Remember that hann window reduces the amplitude of the signal to half in
%frequency domain and that is why scaling is done by N/4 instead of N/2
fft_signal = fft(quantize_signal.*hann(N)',N)/((N/4)*no_quantiz_levels);
freq_scal = linspace (0, 0.5, N/2); % for normalized frequency

figure(2);
plot(freq_scal, abs(fft_signal(1:N/2)), '-x');
hold on
plot(freq_scal(Nw:Nw+2), abs(fft_signal(fin:fin+2)), 'r*');
title('FFT Plot (N = 8192, Hann window is used)');
xlabel('Frequency (Normalized)');

signal_value = (sum((abs(fft_signal(fin:fin+2))).^2))

noise_bins = [fft_signal(3:fin-1) fft_signal(fin+3:N/2)];
% noise_pwr = sum(noise_bins).*conj(noise_bins)))

noise = (sum((abs(noise_bins)).^2))% / ((N/2)-4)

noise_formula = ((1/power(2,adc_resolution))^2)/12

SNR = 10*log10(signal_value/noise)
SNR_byQuantizationFormula = 10*log10(signal_value/noise_formula)


I am anxiously waiting for your reply.

/Iffe

I also noticed that you are only taking half of the signal power and half of the noise power (two errors that cancel each other so, no harm there).

 

[JoannesPaulus]

Below is how I modified your code. As you can see, your "SNR_byQuantizationFormula" gives an incorrect result for a 16 bit quantizer (6.02*adc_resolution+1.76=98.1dB, your formula gives 107.1dB).My result,95.1dB, is closer but not correct yet.

close all;
clear all;

% Generate a simulation signal
fin = 257;
fs = 8192;
N = 8192;
adc_resolution = 16; % Or Quantizer resolution

Nw = floor((fin*N)/fs);
amp = 0.499;
signal = amp*sin(2*pi*Nw/N*[0:N-1]);

figure(1);
subplot (2,1,1);
plot(signal);
title('Original Signal');
ylabel('Signal Amplitude');
xlabel('Time');

% Scalar quantization is implemented over here
no_quantiz_levels = power(2, adc_resolution);
quantize_signal=floor((no_quantiz_levels-1)*signal)/(no_quantiz_levels/2);

subplot (2,1,2);
plot(signal-quantize_signal/2, 'g*');
title('Quantization error');

% FFT and SNR measurement
% Remember that hann window reduces the amplitude of the signal to half in
% frequency domain and that is why scaling is done by N/4 instead of N/2
window_type=ones(1,N);%hann(N)'; %Coherent sampling, no need for windowing
window_scaling_factor=sum(window_type);
fft_signal = fft(quantize_signal.*window_type,N)*(2/window_scaling_factor);
freq_scal = linspace (0, 0.5, N/2); % for normalized frequency
X=fft_signal.*conj(fft_signal);

figure(2);
plot(freq_scal, 10*log10(X(1:N/2)), '-x');
hold on
plot(freq_scal(Nw:Nw+2), 10*log10(abs(fft_signal(Nw:Nw+2))), 'r*');
title('FFT Plot (N = 8192, Hann window is used)');
xlabel('Frequency (Normalized)');

signal_value = (sum(fft_signal(Nw:Nw+2).*conj(fft_signal(Nw:Nw+2))))

noise_bins = [fft_signal(2:fin-1) fft_signal(fin+3:N/2)];
noise = (sum((abs(noise_bins)).^2)) / ((N/2)-3)
noise = 2*(sum(noise_bins.*conj(noise_bins)))
noise_formula = ((1/power(2,adc_resolution))^2)/12

SNR = 10*log10(signal_value/noise)
SNR_byQuantizationFormula = 10*log10(signal_value/noise_formula)

I will keep looking, as I have some spare time...

---------- Post added at 08:57 PM ---------- Previous post was at 08:43 PM ----------

As I copied and pasted the code, I notices that I left a "2*" in the noise calculation:

noise = 2*(sum(noise_bins.*conj(noise_bins)))

That "2*" is redundant because the signal power is calculated on half of the spectrum... so, the SNR is actually correct: 98.1dB!

---------- Post added at 09:23 PM ---------- Previous post was at 08:57 PM ----------

Here is the last bit of information.
Using your notation, the noise formula is incorrect: the rms quantization error is Pnoise_rms=VLSB/sqrt(12), where VLSB is VREF/2^N, VREF=1, in your case.
The rms power of a full-scale sinewave is: Psignal_rms=VREF/(2*sqrt(2)). The SNR is than calculated as 20*log10(Psignal_rms/Pnoise_rms).
In matlab:

noise_formula = (1/2^adc_resolution)/sqrt(12)
SNR_byQuantizationFormula = 20*log10((2*amp/sqrt(8))/noise_formula)

I hope this helps!

 

[iffe]

Hi,

Thanks a lot for such a wonderfull help. Yeh, now I have fixed almost all my bugs by comparing with your code.
Offs......, it seems to be difficult work for me. But I have to do it at all. Can you please recommend me some book or some good references so that I can learn the fundamental concepts for these?

Actually, I have to test a high resolution ADC and I am trying to understand all its parameters with matlab simulations. I am from computer engineering background with some experience in digital.

Also, I have still a problem related to windowing. I have to apply windowing with coherence condition as well. But when I try to do that the SNR starts decreasing with increase in resolution of quantizer even if in your code as well. Upto 12 bits its 4 dB less but in the next it decreases to 85dB for 16 bits which is quite not ok. Can you please let me know that what is the problems in this?

/Iffe

Below is how I modified your code. As you can see, your "SNR_byQuantizationFormula" gives an incorrect result for a 16 bit quantizer (6.02*adc_resolution+1.76=98.1dB, your formula gives 107.1dB).My result,95.1dB, is closer but not correct yet.

close all;
clear all;

% Generate a simulation signal
fin = 257;
fs = 8192;
N = 8192;
adc_resolution = 16; % Or Quantizer resolution

Nw = floor((fin*N)/fs);
amp = 0.499;
signal = amp*sin(2*pi*Nw/N*[0:N-1]);

figure(1);
subplot (2,1,1);
plot(signal);
title('Original Signal');
ylabel('Signal Amplitude');
xlabel('Time');

% Scalar quantization is implemented over here
no_quantiz_levels = power(2, adc_resolution);
quantize_signal=floor((no_quantiz_levels-1)*signal)/(no_quantiz_levels/2);

subplot (2,1,2);
plot(signal-quantize_signal/2, 'g*');
title('Quantization error');

% FFT and SNR measurement
% Remember that hann window reduces the amplitude of the signal to half in
% frequency domain and that is why scaling is done by N/4 instead of N/2
window_type=ones(1,N);%hann(N)'; %Coherent sampling, no need for windowing
window_scaling_factor=sum(window_type);
fft_signal = fft(quantize_signal.*window_type,N)*(2/window_scaling_factor);
freq_scal = linspace (0, 0.5, N/2); % for normalized frequency
X=fft_signal.*conj(fft_signal);

figure(2);
plot(freq_scal, 10*log10(X(1:N/2)), '-x');
hold on
plot(freq_scal(Nw:Nw+2), 10*log10(abs(fft_signal(Nw:Nw+2))), 'r*');
title('FFT Plot (N = 8192, Hann window is used)');
xlabel('Frequency (Normalized)');

signal_value = (sum(fft_signal(Nw:Nw+2).*conj(fft_signal(Nw:Nw+2))))

noise_bins = [fft_signal(2:fin-1) fft_signal(fin+3:N/2)];
noise = (sum((abs(noise_bins)).^2)) / ((N/2)-3)
noise = 2*(sum(noise_bins.*conj(noise_bins)))
noise_formula = ((1/power(2,adc_resolution))^2)/12

SNR = 10*log10(signal_value/noise)
SNR_byQuantizationFormula = 10*log10(signal_value/noise_formula)

I will keep looking, as I have some spare time...

---------- Post added at 08:57 PM ---------- Previous post was at 08:43 PM ----------

As I copied and pasted the code, I notices that I left a "2*" in the noise calculation:

noise = 2*(sum(noise_bins.*conj(noise_bins)))

That "2*" is redundant because the signal power is calculated on half of the spectrum... so, the SNR is actually correct: 98.1dB!

---------- Post added at 09:23 PM ---------- Previous post was at 08:57 PM ----------

Here is the last bit of information.
Using your notation, the noise formula is incorrect: the rms quantization error is Pnoise_rms=VLSB/sqrt(12), where VLSB is VREF/2^N, VREF=1, in your case.
The rms power of a full-scale sinewave is: Psignal_rms=VREF/(2*sqrt(2)). The SNR is than calculated as 20*log10(Psignal_rms/Pnoise_rms).
In matlab:

noise_formula = (1/2^adc_resolution)/sqrt(12)
SNR_byQuantizationFormula = 20*log10((2*amp/sqrt(8))/noise_formula)

I hope this helps!

 

[JoannesPaulus]

Windows are tricky...
The problem you are seeing is power leakage due to the window.
In coherent sampling, the power of the signal is concentrated in one bin.

When you apply a window function, you will notice that the signal power is spread over several bins. Try to plot the DFT of the hann window. You will see that, for a 4096 points DFT, the first point is at -6dB, the second is at -87.8dB, the third at -96.3dB... This is about the power of the quantization noise for a 16-bit ADC!
In order to improve the window performance, you need to increase the number of captured points.

Last, don't forget the 0 frequency! When you calculate the noise power, you will need to skip a few bins near the 0 frequency as well since that is also spread...

Don't forget to push the "helpful" thumb!