How to add variable High Frequency Noise and Baseline drift to signals in a dataset

5 views (last 30 days)
I have a dataset of ecg signals. I have load it all of them into the workspace and filtered with a band-pass butterworth filter. Now I have to add high frequency noise and baseline drift to these filtered signals. The total noise added (high frequency noise + drift) needs to be variable and not be the same.
For the noise can I use a Gaussian or in this case it is wrong? I tried to use awgn, but if it wrong in this context, tell me.
How can I do it?
This is my code:
%% load data and filtering the signals
clear; clc; % clear workspace
Fs = 360; % [Hz]
time_interval = 10; % [sec]
files = dir(fullfile("dataset/","*.mat")); % all dataset files
[b,a] = butter(3,[1 30]/(Fs/2),"bandpass"); % bandpass butterworth filter
filt_sn = zeros(numel(files), 3600); % prealloc
%--load all the data and filter it--%
for i = 1:numel(files)
data(i) = load(fullfile("dataset/",files(i).name)); % load all data
filt_sn(i,:) = filtfilt(b,a,data(i).val); % filter all signals
%--Adding Noise and Baseline Drift--%
% prealloc for noise %
noisy_sn = zeros(numel(files), 3600); % prealloc
SNR = zeros(1,numel(files)); %prealloc
% prealloc for drift
dritf_sn = zeros(numel(files), 3600); % prealloc
A = zeros(1,numel(files)); % prealloc
N = zeros(1,numel(files)); % prealloc
minDriftOffset = zeros(1,numel(files)); % prealloc
maxDriftOffset = zeros(1,numel(files)); % prealloc
for j = 1:numel(files)
%--add noise--%
SNR(j) = 8 + (15 - 8).*rand(1);
noisy_sn(j,:) = awgn(filt_sn(j,:),SNR(j),'measured');
%--add baseline drift--%
%--generate random variable for all signals--%
X = linspace(0,2*pi,length(filt_sn(j,:)));
A(j) = -200 + (200 - (-200)).*rand(1);
N(j) = 50 + (10 - 50).*rand(1);
minDriftOffset(j) = -100 + (0 - (-100)).*rand(1);
maxDriftOffset(j) = 10 + (150 - 10).*rand(1);
%--sum up to distort the signal--%
cos_wave = A(j)*cos(X(j)./N(j));
slantedLine = linspace(minDriftOffset(j),maxDriftOffset(j),length(filt_sn(j,:)));
dritf_sn(j,:) = noisy_sn(j,:)+ cos_wave + slantedLine;
%--plot the signal/s you want--%
%--to plot: chanhe index 'i' in (i,:);
% figure(); subplot(2,1,1); plot(T,noisy_sn(1,:)); grid on;
% title("Noisy ECG signal"); xlabel("time [sec]"); ylabel("voltage [mV]");
% figure(); subplot(2,1,1); plot(T,drfit_sn(1,:)); grid on;
% title("Distorted ECG signal (noise + drift)"); xlabel("time [sec]"); ylabel("voltage [mV]");
  1 Comment
Eamon Gekakis
Eamon Gekakis on 17 Jun 2022
could injecting noise with a sin curve + specified high frequency sample rate + randi() generated noise suffice? Drift can be added by multipling your data by a deviation slope curve.

Sign in to comment.

Accepted Answer

Image Analyst
Image Analyst on 18 Jun 2022
There are lots of kinds of noise. Theoretically what kind does your signal typically have? You can use rand() and randn(). You can also change the amplitude and range of the noise on an element-by-element basis if you want. Do you want to do that, or is your noise relatively stationary across your signal?
For the baseline moving, you have to decide how you want the baseline to move. Do you just want a linear ramp added to the signal? Do you want a sinusoid added to the signal? Do you want the baseline to move fairly randomly, like 1-D Perlin noise was applied to it?
Image Analyst
Image Analyst on 19 Jun 2022
That's basically adding 60 Hz "power line" noise to the fft spectrum (instead of the time domain signal) which doesn't make sense. Why would a pure deterministic signal like that exist in the frequency domain? If you just wanted to add line noise to the signal you'd either add that since wave to the time domain signal, or you'd make a "spike" at the one single element that represented 60 Hz.
But maybe if you want high frequency noise more than 100 Hz, you can try this:
clc; % Clear the command window.
close all; % Close all figures (except those of imtool.)
clear; % Erase all existing variables. Or clearvars if you want.
workspace; % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 15;
markerSize = 30;
fprintf('Beginning to run %s.m ...\n', mfilename);
%--Load Real signal--%
load ('100m.mat');
ecg = val/200;
Fs = 360; % sample frequency [Hz]
L = length(ecg); % Signal Length
T = linspace(0,L/Fs,L); % time axis
F = linspace(-Fs/2, Fs/2, L); % Frequency axis
subplot(3, 2, 1);
plot(T, ecg, 'b-');
grid on;
title('ECG Signal', 'FontSize', fontSize);
xlabel('time (sec)', 'FontSize', fontSize);
ylabel('voltage (mV)', 'FontSize', fontSize);
% Save y axis limits for later.
yl = ylim;
% Get the spectrum of the original signal.
ftEcg = fft(ecg);
shiftedFT = fftshift(abs(ftEcg));
subplot(3, 2, 2);
plot(F, shiftedFT, 'b-');
grid on;
title('Spectrum of original ECG Signal', 'FontSize', fontSize);
xlabel('Frequency (Hz)', 'FontSize', fontSize);
ylabel('Signal', 'FontSize', fontSize);
% Make noise vector
amplitude = 2;
noiseSpectrum = amplitude * rand(1, length(T));
shiftedNS = fftshift(noiseSpectrum);
% Zero out less than 100 Hz.
index = round(length(shiftedNS) * 100 / Fs)
middleIndex = length(shiftedNS) / 2
shiftedNS((middleIndex - index) : (middleIndex + index)) = 0;
subplot(3, 2, 3);
plot(F, shiftedNS, 'b-');
grid on;
title('Spectrum of Noise Signal', 'FontSize', fontSize);
xlabel('Frequency (Hz)', 'FontSize', fontSize);
ylabel('Signal', 'FontSize', fontSize);
% Add noise spectrum to original signal spectrum
noisySpectrum = abs(shiftedFT + shiftedNS);
subplot(3, 2, 4);
plot(T, noisySpectrum, 'b-');
grid on;
title('Noisy ECG Spectrum', 'FontSize', fontSize);
xlabel('Frequency (Hz)', 'FontSize', fontSize);
ylabel('signal', 'FontSize', fontSize);
% Shift it back to get the spectrum in a form where we can use ifft:
noisySpectrum = ifftshift(noisySpectrum);
subplot(3, 2, 5);
plot(T, noisySpectrum, 'b-');
grid on;
title('Noisy ECG Spectrum', 'FontSize', fontSize);
xlabel('Frequency (Hz)', 'FontSize', fontSize);
ylabel('signal', 'FontSize', fontSize);
% Return to time domain
noisyTimeDomainSignal = ifft(noisySpectrum);
subplot(3, 2, 5:6);
plot(T, noisyTimeDomainSignal, 'b-');
grid on;
title('Noisy ECG Signal', 'FontSize', fontSize);
xlabel('time (sec)', 'FontSize', fontSize);
ylabel('voltage (mV)', 'FontSize', fontSize);
No guarantees - you'll have to review my math.

Sign in to comment.

More Answers (0)

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!