MATAB Code
%% ==========================================================
% WSS / Noise Demonstration Simulator
% Sinusoid + Multiple Noise Types
% Distribution + Autocorrelation + PSD
% LTI System Demonstration
%% ==========================================================
clear;
close all;
clc;
%% Parameters
Fs = 1000; % Sampling Frequency
N = 5000; % Number of Samples
t = (0:N-1)/Fs;
A = 1; % Sinusoid Amplitude
f0 = 20; % Sinusoid Frequency
signal = A*sin(2*pi*f0*t);
noiseVariance = 0.25;
noiseStd = sqrt(noiseVariance);
%% ==========================================================
% Generate Noise Types
%% ==========================================================
% Gaussian Noise
gaussianNoise = noiseStd*randn(1,N);
% Uniform Noise (same variance)
uniformNoise = (rand(1,N)-0.5)*sqrt(12*noiseVariance);
% Laplace Noise
u = rand(1,N)-0.5;
b = sqrt(noiseVariance/2);
laplaceNoise = -b*sign(u).*log(1-2*abs(u));
% Binary Noise
binaryNoise = noiseStd*sign(rand(1,N)-0.5);
% Pink Noise (simple approximation)
whiteNoise = randn(1,N);
pinkNoise = filter(1,[1 -0.98],whiteNoise);
pinkNoise = pinkNoise/std(pinkNoise)*noiseStd;
%% ==========================================================
% Create Noisy Signals
%% ==========================================================
sigGaussian = signal + gaussianNoise;
sigUniform = signal + uniformNoise;
sigLaplace = signal + laplaceNoise;
sigBinary = signal + binaryNoise;
sigPink = signal + pinkNoise;
%% ==========================================================
% Figure 1 : Time Domain Signals
%% ==========================================================
figure('Name','Sinusoid with Different Noise Types',...
'Position',[100 100 1200 800]);
subplot(3,2,1)
plot(t(1:1000),signal(1:1000),'k','LineWidth',1.5)
title('Pure Sinusoid')
xlabel('Time (s)')
ylabel('Amplitude')
grid on
subplot(3,2,2)
plot(t(1:1000),sigGaussian(1:1000))
title('Sinusoid + Gaussian Noise')
grid on
subplot(3,2,3)
plot(t(1:1000),sigUniform(1:1000))
title('Sinusoid + Uniform Noise')
grid on
subplot(3,2,4)
plot(t(1:1000),sigLaplace(1:1000))
title('Sinusoid + Laplace Noise')
grid on
subplot(3,2,5)
plot(t(1:1000),sigBinary(1:1000))
title('Sinusoid + Binary Noise')
grid on
subplot(3,2,6)
plot(t(1:1000),sigPink(1:1000))
title('Sinusoid + Pink Noise')
grid on
%% ==========================================================
% Figure 2 : Noise Histograms
%% ==========================================================
figure('Name','Noise Distributions',...
'Position',[150 100 1200 800]);
subplot(3,2,1)
histogram(gaussianNoise,50,'Normalization','pdf')
title('Gaussian Distribution')
grid on
subplot(3,2,2)
histogram(uniformNoise,50,'Normalization','pdf')
title('Uniform Distribution')
grid on
subplot(3,2,3)
histogram(laplaceNoise,50,'Normalization','pdf')
title('Laplace Distribution')
grid on
subplot(3,2,4)
histogram(binaryNoise,50,'Normalization','pdf')
title('Binary Distribution')
grid on
subplot(3,2,5)
histogram(pinkNoise,50,'Normalization','pdf')
title('Pink Noise Distribution')
grid on
%% ==========================================================
% Figure 3 : Distribution Comparison
%% ==========================================================
figure('Name','Distribution Comparison');
hold on
[f,x] = ksdensity(gaussianNoise);
plot(x,f,'LineWidth',2)
[f,x] = ksdensity(uniformNoise);
plot(x,f,'LineWidth',2)
[f,x] = ksdensity(laplaceNoise);
plot(x,f,'LineWidth',2)
legend('Gaussian','Uniform','Laplace')
xlabel('Amplitude')
ylabel('Probability Density')
title('Noise Probability Density Functions')
grid on
%% ==========================================================
% Figure 4 : Autocorrelation
%% ==========================================================
figure('Name','Noise Autocorrelation',...
'Position',[150 100 1200 800]);
subplot(3,2,1)
autocorr(gaussianNoise,100)
title('Gaussian')
subplot(3,2,2)
autocorr(uniformNoise,100)
title('Uniform')
subplot(3,2,3)
autocorr(laplaceNoise,100)
title('Laplace')
subplot(3,2,4)
autocorr(binaryNoise,100)
title('Binary')
subplot(3,2,5)
autocorr(pinkNoise,100)
title('Pink')
%% ==========================================================
% Figure 5 : Power Spectral Density
%% ==========================================================
figure('Name','Power Spectral Density',...
'Position',[150 100 1200 800]);
subplot(3,2,1)
pwelch(gaussianNoise)
title('Gaussian PSD')
subplot(3,2,2)
pwelch(uniformNoise)
title('Uniform PSD')
subplot(3,2,3)
pwelch(laplaceNoise)
title('Laplace PSD')
subplot(3,2,4)
pwelch(binaryNoise)
title('Binary PSD')
subplot(3,2,5)
pwelch(pinkNoise)
title('Pink PSD')
%% ==========================================================
% Figure 6 : WSS Through LTI System
%% ==========================================================
h = [0.2 0.5 0.2];
x = gaussianNoise;
y = conv(x,h,'same');
figure('Name','WSS Through LTI');
subplot(2,2,1)
plot(x(1:500))
title('Input WSS Process')
grid on
subplot(2,2,2)
plot(y(1:500))
title('Output Process')
grid on
subplot(2,2,3)
[Rx,lags] = xcorr(x,100,'biased');
plot(lags,Rx,'LineWidth',1.5)
title('Input Autocorrelation')
xlabel('Lag')
grid on
subplot(2,2,4)
[Ry,lags] = xcorr(y,100,'biased');
plot(lags,Ry,'LineWidth',1.5)
title('Output Autocorrelation')
xlabel('Lag')
grid on
%% ==========================================================
% Figure 7 : WSS vs Non-WSS
%% ==========================================================
xWSS = randn(1,N);
xNonWSS = (1 + 0.002*(1:N)).*randn(1,N);
window = 200;
runningMeanWSS = movmean(xWSS,window);
runningMeanNonWSS = movmean(xNonWSS,window);
figure('Name','WSS vs Non-WSS');
subplot(2,1,1)
plot(runningMeanWSS,'LineWidth',1.5)
title('Running Mean of WSS Process')
ylabel('Mean')
grid on
subplot(2,1,2)
plot(runningMeanNonWSS,'r','LineWidth',1.5)
title('Running Mean of Non-WSS Process')
ylabel('Mean')
xlabel('Sample Index')
grid on
%% ==========================================================
% Figure 8 : Correlation Matrix Visualization
%% ==========================================================
[R,lags] = xcorr(x,50,'biased');
Rpositive = R(lags>=0);
figure('Name','Correlation Matrix');
imagesc(toeplitz(Rpositive))
colormap(jet)
colorbar
title('Autocorrelation Matrix of WSS Process')
xlabel('Sample Index')
ylabel('Sample Index')
%% ==========================================================
% Summary
%% ==========================================================
disp('========================================');
disp('Simulation Complete');
disp('Generated:');
disp('1. Noisy sinusoid signals');
disp('2. Noise distributions');
disp('3. PDF comparison');
disp('4. Autocorrelation plots');
disp('5. Power spectral densities');
disp('6. WSS -> LTI -> WSS demonstration');
disp('7. WSS vs Non-WSS comparison');
disp('8. Correlation matrix visualization');
disp('========================================');
Output
