In this post, i am going to show how to add noise to the time series signal (assume the time series signal is sech(t) in this post) via its frequency domain spectrum in MATLAB
The MATLAB code below plots the time domain (i.e., u) and frequency domain (i.e., ut) of signal sech(t) (using 512 points in time domain data with periodic domain set to T=30)
clear all; close all; clc; T=30; n=512; t2=linspace(-T/2, T/2, n+1); t=t2(1:n); k=(2*pi/T)*[0:n/2-1 -n/2:-1]; ks=fftshift(k); u=sech(t); ut=fft(u); subplot(2, 1, 1); plot(t, u); subplot(2, 1, 2); plot(ks, abs(fftshift(ut))/max(abs(fftshift(ut)))); axis([-25 25 0 1]);
The code abs(fftshift(ut))/max(abs(fftshift(ut))) normalize the power spectrum value to be between 0 and 1
The code below shows the resulting plot
The code below adds a white noise with normal distribution (noise~N(0, 1)) to the frequency domain (i.e., utn) which is then converted to back to the time domain (i.e. un)
clear all; close all; clc; T=30; n=512; t2=linspace(-T/2, T/2, n+1); t=t2(1:n); k=(2*pi/T)*[0:n/2-1 -n/2:-1]; ks=fftshift(k); u=sech(t); ut=fft(u); noise=1; utn=ut+noise*(randn(1, n)+i*randn(1, n)); un=ifft(utn); subplot(2, 1, 1); plot(t, u, 'k', t, abs(un), 'm'); subplot(2, 1, 2); plot(ks, abs(fftshift(ut))/max(abs(fftshift(ut))), 'k', ks, abs(fftshift(utn))/max(abs(fftshift(utn))), 'm'); axis([-25 25 0 1]);
The image below shows the resulting plot
No comments:
Post a Comment