Calculate the FFT of an ECG signal

82 次查看(过去 30 天)
Hi all,
currently i'am trying to transform an ecg signal into frequency domain.
I use an ecg signal from MIT-BIH Arrhythmia Database (physionet). To read the ecg file i use the .dat file reader from Said Bourezg (https://de.mathworks.com/matlabcentral/fileexchange/49822-open_ecg-ecg-dat-file-reader).
After reading the signal i add a noisy sine wave to corrupt the signal. Finally I want to calculate the Fourier transform of the noisy signal.
The ECG signal is sampled at 150 Hz and has a duration of 10 seconds (3500 samples/sec).
I would like to know if my frequency plot is correct. It looks a bit weird. It has a relatively high peak at 1 Hz and no other (or small ) magnitudes. And of course you can see the 50 Hz noise. I was actually expecting more frequencies.
clear all
clc
close all
[filename, pathname] = uigetfile('*.dat', 'Open file .dat');
if isequal(filename, 0) || isequal(pathname, 0)
disp('File input canceled.');
ECG_Data = [];
else
fid=fopen(filename,'r');
end;
time=10; % duration of ecg signal
fs = 150; %sample frequency
f=fread(fid,2*fs*time,'ubit12');
len = length(f(1:2:length(f))); %length of the signal
time_step = 1/fs;
max_time = len/fs; % duration of your signal in seconds depending on sampling rate
t = time_step:time_step:max_time; % this is our time vector.
Orig_Sig=f(1:2:length(f))/len; %compute the clean ecg signal
subplot(4,1,1);
plot(t,Orig_Sig)
noise_coeff = 0.1; %randomize the amplitude of the noisy ecg signal
noise_signal = 5*sin(2*pi*50*t); %add 50Hz noise from power supply
noise_signal = noise_signal'; %transpose vector
dirty_signal = Orig_Sig + noise_coeff*noise_signal; %mix the original ecg signal with the noisy 50 Hz sine
subplot(4,1,2);
plot(t,noise_signal);
subplot(4,1,3);
plot(t,dirty_signal);
NFFT = 2 ^ nextpow2(length(dirty_signal)); %compute FFT length depends on the signal length
Y = fft(dirty_signal,NFFT); %compute the fft of the noisy signal
Y = Y(1:NFFT/2); %we only need a one sided fft plot
Y_abs = 1/NFFT*abs(Y); %calculate the magnitude and normalize the spectrum
f_fft = (0:NFFT/2-1)*fs/NFFT; %scale the frequency axe and calculate the coresponding frequencys
subplot(4,1,4);
plot(f_fft,Y_abs);
I also added the used ecg signal (.dat file) and a frequency plot.
If anyone finds a mistake I would be very grateful.
Best regards.

采纳的回答

Daniel M
Daniel M 2019-10-8
编辑:Daniel M 2019-10-8
I haven't ran your code but the plot looks pretty normal. The peak you see around 0 Hz is a DC component added to your signal. You can subtract the mean of your data using detrend(), and see if that helps. Otherwise, try looking at the plot using semilogy. The resolution would be better if you had more data than just 10 seconds.
  3 个评论
Nur Shahirah
Nur Shahirah 2020-4-14
编辑:Nur Shahirah 2020-5-9
Can i know in which line you put the this code? i try to figure it out, yet i still have problems with the frequency graph. Thanks in advance.
9/5/2020 : I still can't configure this. Help?

请先登录,再进行评论。

更多回答(0 个)

Community Treasure Hunt

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

Start Hunting!

Translated by