Help filtering a signal using fft
50 views (last 30 days)
Show older comments
Hey,
I'm trying to filter a signal with the use of the fft. I would like to keep the lower frequencies and cut the higher frequencies.
I have given the code below.
I assume that I'm wrong in defining my filter variable... Right now, my FFT has length ~8000, and I wanted to keep the first 300 samples of that fft, and cut off the rest.
clear all;
close all;
load('dataTP.mat');
L=length(dataTP.data);
%create data vector, the one I want to filter
for i=1:length(dataTP.data)
data(i)=sqrt((dataTP.data(i,1))^2+(dataTP.data(i,2))^2+(dataTP.data(i,3))^2);
end
data = data-mean(data);
%create filter for use in fftfilt
filter=zeros(1,length(data));
filter(1:300)=1;
filtered_data=fftfilt(filter, data);
%create fft to visualize in graph
data_fft=fft(data);
P2 = abs(data_fft/L);
P1 = P2(1:L/2+1);
P1(2:end-1) = 2*P1(2:end-1);
f = 200*(0:(L/2))/L;
%plot fft visualization
figure(1)
plot(f,P1)
%plot raw data
figure(2)
plot(data);
%plot filtered data vs raw data
figure(3)
plot(F)
hold on
plot(data, 'r')
hold off
My first attempt was made like this:
data_fft=fft(data);
filter=zeros(1,length(data));
filter(1:300)=1;
filtered_fft=data_fft.*filter;
filtered_data=ifft(filtered_fft);
but that didn't work too well either... Hope you can help!
I have attached the data File.
Thanks!
2 Comments
John BG
on 9 Oct 2016
could please make available the dataTP1.mat file, or a portion of it so the readers of your question can reproduce your set-up?
Accepted Answer
Star Strider
on 9 Oct 2016
Ideally, you should also have a time vector. This allows you to define specific frequencies rather than referring to elements of the Fourier transform vector. I created a time vector, so you can substitute it for yours or leave my code as it is.
See if this does what you want:
d = load('simon DataTP.mat');
DataRaw = d.dataTP.data;
data = sqrt(sum(DataRaw.^2,2)); % Euclidean
Ts = 1; % Sampling Interval
Fs = 1/Ts; % Sampling Frequency
Fn = Fs/2; % Nyquist Frequency
L = size(data,1); % Length Of ‘data’ Vector
t = 1:L*Ts; % Time Vector
FTdata = fft(data)/L; % Fourier Transform
Fv = linspace(0, 1, fix(L/2)+1)*Fn; % Frequency Vector (One-Sided FFT)
Iv = 1:length (Fv); % Index Vector
Iv300 = 1:300; % ‘First 300’ FFT Frequencies
F300 = Fv(300); % Frequency Corresponding To 300th Element
figure(1)
plot(Fv, abs(FTdata(Iv))*2)
grid
figure(2)
plot(Fv(Iv300), abs(FTdata(Iv300))*2)
grid
FLen = 48; % Discrete Filter Order
b_filt = fir1(FLen, F300/Fn); % Design FIR Filter
figure(3)
freqz(b_filt, 1, 4096, Fs)
data_filtered = fftfilt(b_filt, data);
figure(4)
subplot(2,1,1)
plot(t, data)
title('Original Data')
grid
subplot(2,1,2)
plot(t, data_filtered)
title('Filtered Data')
grid
2 Comments
Star Strider
on 10 Oct 2016
My pleasure!
Overnight, my primary HP Win 8.1 MATLAB computer died (it needs yet another power management chip, as it seems to need every couple years), and this Win 10 computer (running R2015b) won't let me load .mat files with either Firefox or IE. (My opinion of Win 10 is not fit for an open forum.) So I can’t load your file to test it. I'll access Answers on another machine or on my phone to see if I can download it and then transfer it via Bluetooth later. I’m tired of having to do work-arounds for Win 10’s myriad deficiencies.
See if this works for ‘t’:
t = linspace(0, 1, L)*Ts;
More Answers (0)
See Also
Categories
Find more on Digital Filter Analysis in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!