Issue with FFT/Power Spectral Density

5 次查看(过去 30 天)
Can someone explain what I'm doing wrong here? I wanna evaluate the PSD of ECG data (electrocardiogram), but I end up with this enormous peak in the fourier domain that I don't know what to do with. Where does this come from? What to do about it?
clear
load ecg_stress_test.mat
n=length(ecg_data);
fhat=fft(ecg_data,n);
PSD=fhat.*conj(fhat)/n;
figure
subplot(2,1,1), plot(ecg_data), xlim([0 length(ecg_data)])
subplot(2,1,2), plot(PSD), xlim([0 length(PSD)])

采纳的回答

Star Strider
Star Strider 2021-11-5
Subtract the mean of the EKG data before calculating the fft. Also, it likely needs to be rescaled, because the amplitude should be between ±1 mV.
The signal needs to be filtered. Some of the noise is broadband, however much can be filtered with a frequency-selective filte —
LD = load('Sebastian Daneli ecg_stress_test.mat')
LD = struct with fields:
ecg_data: [450001×1 double] fs: 250
EKG = LD.ecg_data;
Fs = LD.fs;
L = numel(EKG);
t = linspace(0, 1, L)/Fs;
[EKGfilt,df] = bandpass(EKG, [1 25], Fs, 'ImpulseResponse','iir');
figure
subplot(3,1,1)
plot(t, EKG)
grid
subplot(3,1,2)
idx = 200000:201000;
plot(t(idx), EKG(idx))
grid
subplot(3,1,3)
idx = 200000:201000;
plot(t(idx), EKGfilt(idx))
grid
Fn = Fs/2;
NFFT = 2^nextpow2(L);
FTEKGfilt = fft(EKGfilt,NFFT)/L;
Fv = linspace(0, 1, NFFT/2+1)*Fn;
Iv = 1:numel(Fv);
[ampmax,idx] = max(abs(FTEKGfilt(Iv))*2);
fprintf('Amplitude = %.2f mV at %.2f Hz',ampmax,Fv(idx))
Amplitude = 1.20 mV at 1.92 Hz
figure
plot(Fv, abs(FTEKGfilt(Iv))*2)
grid
The pwelch function is probably the mose frequently used method to determine the power spectral density. To see the the signal spectrum and power spectrum as a function of time, use the pspectrum or spectrogram functions.
.

更多回答(1 个)

Paul
Paul 2021-11-4
That peak is at dc because the mean of the data is very large. Try
fhat=fft(detrend(ecg_data),n);
to take out the mean before taking the fft.
  6 个评论
Paul
Paul 2021-11-8
I was confused abou the use of detrend. It should have been
fhat = fhat=fft(detrend(ecg_data,'constant'),n);
to subtract the mean.

请先登录,再进行评论。

类别

Help CenterFile Exchange 中查找有关 Spectral Measurements 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by