Unexpected power spectrum of a signal with uniform random amplitude values
13 views (last 30 days)
Show older comments
Vlad Atanasiu
on 5 Jan 2022
Commented: William Rose
on 11 Jan 2022
What should the power spectrum of a signal with uniform random amplitude values be? I expect it to be uniform within the constraints of a finite signal, with no frequency dominating.
However, using the following code I do not obtain the expected quasi-uniform power spectrum. In the attached diagram file, the magnitude ressembles a gamma function. Why is this so?
![Spectrum of signal with uniform random amplitude values](https://www.mathworks.com/matlabcentral/answers/uploaded_files/854340/Spectrum%20of%20signal%20with%20uniform%20random%20amplitude%20values.png)
% Spectrum of signal with uniform random amplitude values
clear all
% signal
n = 2^9; % 2^9 = 512, 2^13 = 8192 % signal length
t = 1:n; % time
y = unifrnd(-1,1,[1,n]); % amplitude, uniform random
% mu = 0;
% sigma = 1;
% y = random('Normal',mu,sigma,[1,n]); % amplitude, normal random
% spectral analysis
nyq = ceil(n/2); % Nyquist frequency
f = fft(y); % Fourier transform
f = f(1:nyq); % half the symmetrical spectrum
m = abs(f); % magnitude
m = m.^2; % power
%m = log(m.^2); % log power
%m = sqrt(m);
p = angle(f); % phase
% sorted values
ys = sort(y,'descend');
ms = sort(m,'descend');
ps = sort(p,'descend');
% histograms
hn = 100;
hx = 1:hn;
hy = histcounts(y,hn);
hm = histcounts(m,hn);
hp = histcounts(p,hn);
% visualization
figure('Name',['Analysis of noise spectra (n = ',num2str(n),')'])
% signal
subplot(3,2,1)
hold on
plot(t,y,'b')
plot(t,ys,'r','LineWidth',2)
hold off
xlim('padded')
ylim('padded')
axis square
box on
title('Signal')
xlabel({'Time','Red line: sorted amplitude values.'})
ylabel('Amplitude')
subplot(3,2,2)
stairs(hy,hx,'b')
xlim('padded')
ylim('padded')
axis square
title('Histogram')
xlabel('Counts')
ylabel('Amplitude Values')
% magnitude/power
subplot(3,2,3)
hold on
plot(t(1:nyq),m,'b')
plot(t(1:nyq),ms,'r','LineWidth',2)
hold off
%set(gca,'XScale','log')
%set(gca,'YScale','log')
xlim('padded')
ylim('padded')
axis square
box on
%title('Magnitude')
title('Power (Magnitude^2)')
%title('Log Power (Magnitude^2)')
xlabel('Frequency')
%xlabel('Log Frequency')
ylabel('Amplitude')
%ylabel('Log Amplitude')
subplot(3,2,4)
stairs(hm,hx,'b')
xlim('padded')
ylim('padded')
axis square
title('Histogram')
xlabel('Counts')
ylabel('Amplitude Values')
% phase
subplot(3,2,5)
hold on
plot(t(1:nyq),p,'b')
plot(t(1:nyq),ps,'r','LineWidth',2)
hold off
xlim('padded')
ylim('padded')
axis square
box on
title('Phase')
xlabel('Frequency')
ylabel('Radians')
subplot(3,2,6)
stairs(hp,hx,'b')
xlim('padded')
ylim('padded')
axis square
title('Histogram')
xlabel('Counts')
ylabel('Amplitude Values')
0 Comments
Accepted Answer
William Rose
on 10 Jan 2022
You are correct that the power spectrum of a random signal is expected to be approximately flat. However, the magntidues of the spectrum, computed as you have done, with no spectrum averaging, have a standard deviation as large as the mean value. This is explained in textbooks on random signals.
If you do spectrum averaging, you reduce the standard deviation, although you lose spectral resolution. You can demonstrate this using the built-in function pwelch().
You have added a plot of the sorted spectral data to the plot of the magnitude spectrum, at subplot(3,2,3). This red line plot is incorrect, because the y values are plotted versus specific frequencies. But when you sort the magntiudes, the frequencies get all mixed up. So the red plot should not be there.
Power spectral estimates from Gaussian noise have a chi-squared distribution. If you use Gaussian noise instead of uniform, and make a histogram of the spectrum, you will confirm that this is true. I have done it in the past.
Good luck with your work.
4 Comments
More Answers (0)
See Also
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!