Remove peak at 0 hz of fft

101 ビュー (過去 30 日間)
Emmy
Emmy 2020 年 7 月 6 日
コメント済み: Image Analyst 2020 年 7 月 10 日
Hi everybody,
I want to calculate the fft of my data, but I get a large peak at 0 Hz. I have read all the answers to the previous questions. Mostly recommended using x=x-mean(x) or x=detrend(x). I have tried that but nothing works.
Another suggestion was to look at the diff(x) and diff(diff(x)). I tried that and both give a noisy straight line.
I dont now what else to do so can somebody help me?
Thanks in advance!
x=load('data.txt');
% diffx=diff(x);
% diffdiffx=diff(diff(x));
xm=x-mean(x);
xd=detrend(x);
y=detrend(xm,'constant');
N=length(x);
fs=1000;
Ts=1/fs;
f=(0:1/N:1/2-1/N)*fs;
% figure
% plot(x)
% hold on
% plot(diffx)
% plot(diffdiffx)
% hold off
% legend('x','diffx','diffdiffx')
xhat=fft(x,N);
xhat(1)=0;
figure
plot(f,abs(xhat(1:N/2)));
title('Normal data')
xhat=fft(xm,N);
xhat(1)=0;
figure
plot(f,abs(xhat(1:N/2)));
title('Removed mean data')
xhat=fft(xd,N);
xhat(1)=0;
figure
plot(f,abs(xhat(1:N/2)));
title('Detrend data')
xhat=fft(y,N);
xhat(1)=0;
figure
plot(f,abs(xhat(1:N/2)));
title('Removed mean and detrend data')
[S,fm]=periodogram(x,[],N,fs,'power');
figure
plot(fm,10*log10(S))
title('Periodogram')
[S,fm]=periodogram(xm,[],N,fs,'power');
figure
plot(fm,10*log10(S))
title('Periodogram removed mean')
[S,fm]=periodogram(xd,[],N,fs,'power');
figure
plot(fm,10*log10(S))
title('Periodogram detrend data')
[S,fm]=periodogram(y,[],N,fs,'power');
figure
plot(fm,10*log10(S))
title('Periodogram removed mean and detrend data')
[Sw,fw]=pwelch(x,ones(N,1),0,N,fs,'power');
figure
plot(fw,Sw)
title('Welch method')
[Sw,fw]=pwelch(xm,ones(N,1),0,N,fs,'power');
figure
plot(fw,Sw)
title('Welch method removed mean')
[Sw,fw]=pwelch(xd,ones(N,1),0,N,fs,'power');
figure
plot(fw,Sw)
title('Welch method detrend data')
[Sw,fw]=pwelch(y,ones(N,1),0,N,fs,'power');
figure
plot(fw,Sw)
title('Welch method removed mean and detrend data')
  2 件のコメント
Mehmed Saad
Mehmed Saad 2020 年 7 月 6 日
apply a hpf
Emmy
Emmy 2020 年 7 月 6 日
Applying a highpass filter did not solve the problem. There still is a peak from the point of the high pass filter. And it doesn't matter at which point I take the high pass filter. Even when I look in the normal data, where the peak ends, there still is a peak after the high pass filter.

サインインしてコメントする。

採用された回答

David Goodmanson
David Goodmanson 2020 年 7 月 9 日
Hi Emmy
It's all in how you view things. Literally. Here is what happens when using x - mean(x).
fs = 1000;
N = 30000;
figure(1)
plot(x)
xm = x - mean(x);
y = fftshift(abs(fft(xm))); % shift zero frequency to the center
f = (-N/2:N/2-1)*fs/N; % make f array match
y(15001) = 1; % adjust 0 Hz value in fft
figure(2)
semilogy(f,y)
figure(3)
semilogy(f,y)
xlim([-50 50])
Both x and y have 30000 elements. Ordinarily, y(1) would be the element corresponding to zero frequency. For plotting purposes, fftshift puts that element at the center of the array, element 15001.
Since the y values have a large range, it makes sense to look at abs(y) on a semilog plot. Note that since y is the fft of a real function, abs(y) is symmetric about f==0.
After subtracting the mean, the DC value should be zero. Because of 'numerical error' it is small but nonzero. For plot purposes I arbitrarily set that value to 1 so as to not have a large negative spike in the semilog plot.
In figure 2, you can see some lower frequency stuff rising above the noise. In the zoomed-in figure 3. there is something interesting happening between about 17-20 Hz. If you zoom in a little bit in figure 3 you will see no peak at f=0, with the largest values at one freq point on either side, falling off from there.
  1 件のコメント
Image Analyst
Image Analyst 2020 年 7 月 10 日
Adding plots for visualization purposes:
x=load('data.txt');
fs = 1000;
N = 30000;
subplot(3, 1, 1);
plot(x)
grid on;
title('x', 'FontSize', 20);
xm = x - mean(x);
y = fftshift(abs(fft(xm))); % shift zero frequency to the center
f = (-N/2:N/2-1)*fs/N; % make f array match
y(15001) = 1; % adjust 0 Hz value in fft
subplot(3, 1, 2);
semilogy(f,y)
grid on;
title('FT', 'FontSize', 20);
subplot(3, 1, 3);
semilogy(f,y)
xlim([-50 50])
grid on;
title('Zoomed FT', 'FontSize', 20);

サインインしてコメントする。

その他の回答 (1 件)

Mehmed Saad
Mehmed Saad 2020 年 7 月 6 日

FFT of X

FFT of X High pass filtered

カテゴリ

Help Center および File ExchangeDescriptive Statistics についてさらに検索

タグ

Community Treasure Hunt

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

Start Hunting!

Translated by