Determining the ideal smoothing factor

7 views (last 30 days)
David
David on 7 Mar 2016
Commented: Star Strider on 30 Oct 2016
I am researcher using matlab to analise a very large sample of accelerometer data collected using the MPU9150 IC. My first aim is to plot the FFT for this sample.
As I am new to MATLAB I have taken over month to plot the FFT graph but I believe its highly dependent on the smoothing factor. I wish to first illustrate my approach and then showcase my full code for your kind reference. Sampling frequency of of the accelerometer is 10 Hz.
Approach in exact order.
  1. Import the raw acceleration vector with x,y,z axes
  2. Determine butterworth Filter parameters
  3. Apply butterworth filter using filtfilt ()
  4. use smooth () function with a DC smooth factor. ( smoothing factor effects FFT graph significantly)
  5. obtain FFT using abs(fft(acceleration_DC_shifted_Filtered_Vector))
  6. smooth the FFT using Smooth function (this smoothing has no effect on FFT)#plot FFT
%importing data from workspace.
accelerationd2 = Xaccelday2;
% Setting basic Parameters. Number of samples should be same in all day files. number_of_samples = length (accelerationd2); sampling_frequency = 10;
%Setting filtering paramters.
cut_off_frequency = 5;
wn = cut_off_frequency/sampling_frequency;
order_of_butterworth_filter = 1;
[b1,a1]= butter(order_of_butterworth_filter,wn,'low');
%Filtering data from all ten days.
accelerationd2_filtered = filtfilt(b1,a1,accelerationd2);
%DC shift removal
accelerationd2_smoothned= smooth(accelerationd2_filtered,36.5);%acceleration 1filtered was passed into all smoothfunctions
accelerationd2_DC_shifted = accelerationd2_filtered-accelerationd2_smoothned;
%Computing FFT
accelerationd2_FFT = smooth(abs(fft(accelerationd2_DC_shifted)),8000); %1000 is ballpark
bin_vals = (0 : number_of_samples-1);
fax_Hz = bin_vals*sampling_frequency/number_of_samples;
N_2 = ceil(number_of_samples/2);
%Detecting peaks
figure;
plot((fax_Hz(1:N_2)),(accelerationd2_FFT(1:N_2)),'r');
hold on
legend('Day-02');
%[maxtabDay1, mintabDay1] = peakdet(accelerationd1_FFT(1:N_2), 15, fax_Hz(1:N_2));
[maxtabDay2, mintabDay2] = peakdet(accelerationd2_FFT(1:N_2), 15, fax_Hz(1:N_2));
% This is where peaks are dressed with markers. Markers must have the same color of the line. Markers can be different symbols. Types of markers that can be used : ^ , *
plot(maxtabDay2(:,1), maxtabDay2(:,2), 'r*');
xlabel('Frequency (Hz)');
ylabel('Magnitude (dB)');
title('Single-sided Magnitude Spectrum (Hz)');
%hold on
hold off
%axis tight
My problem is that the FFT changes highly as the smoothing factor in
accelerationd2_smoothned= smooth(accelerationd2_filtered,36.5)
is changed. Any help would be so much appreciated.
This is my very first post on this forum. I apologize if I have overlooked any rules.
  2 Comments
Rena Berman
Rena Berman on 30 Oct 2016
(Answers dev) restored question
Star Strider
Star Strider on 30 Oct 2016
Thank you, Rena!
It would be nice for this to be Standard Policy.

Sign in to comment.

Answers (1)

Star Strider
Star Strider on 9 Mar 2016
First, see the R2015a documentation for fft, specifically the code between the first two plot figures. That is everything you need to know about calculating it and plotting it.
Second, you need to re-design your filter. My filter design procecure is in How to design a lowpass filter for ocean wave data in Matlab? Also see my recent bandpass filter design in: How to Butterworth Filter with Bandpass [10 500] with sampling rate 1000?. Note that your sampling frequency is the inverse of the sampling interval:
time = ... time vector ...;
Ts = mean(diff(time));
Fs = 1/Ts;
  2 Comments
David
David on 16 Mar 2016
thank you for your time Star Strider I appreciate your time.
I too have applied butterworth filter. And I see nothing wrong with it. Could you be kind enough to tell me if I have done anything wrong in computing the FFT?
My FFT seems to take a shift in frequency when I alter the smoothing factor. That is my problem. The answer you directed me is about butterworth filtering. I have done the filtering accordingly. My problem is FFT does not work when the acceleration signal is not about the time axis. i.e. when there is a DC shift due to gravity on the acceleration signal. So smooth has been used to remove it.
Star Strider
Star Strider on 16 Mar 2016
My pleasure.
I always recommend using this fft procedure. Not particularly the code between the top two figure plots. It is everything you need to know with respect to calculating the fft of your signal.
I do not have the Curve Fitting Toolbox, so I do not have the smooth function. When I looked it up in the online documentation, I discovered that it is a moving average filter. I discourage the use of moving average filters because they are not designed with any sort of frequency selection. If they are used with the filter function, they will also have the phase distortion inherent in using it instead of the Signal Processing Toolbox filtfilt function, which has no phase distortion. I would not use any sort of smoothing of the Fourier transformed data in any event. That obscures the exact information you want to get from the fft.
If you want to eliminate gravity (a d-c offset), use a band-pass filter with a low-frequency cutoff (stop-band) frequency in the range of 1 to 10 Hz, depending on your sampling frequency. (Higher sampling frequencies will allow lower cutoff frequencies and steeper roll-offs.) You will need the Signal Processing Toolbox for this, but it is well worth having if you do any sort of sophisticated filter design and implementation.

Sign in to comment.

Products

Community Treasure Hunt

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

Start Hunting!