How to filter noise from time-frequency data and find natural frequency of a cantilever?

Hi,

I am studying the natural frequency of a cantilever. For my experiments I collect the deflection data (X and Y coordinates) of a mirror at the tip of the cantilever using a photodetector. I am using Fast Fourier Transform to convert the time domain data into frequency domain. I have attached two sample datasets. The code that I am using for the analysis is:

data=xlsread('sample_data_1.xlsx');

t=data(:,1);

y=data(:,3);

Fs=500;

nfft=2000;

Y=fft(y,nfft);

Y=Y(1:nfft/2);

mx=abs(Y);

f=(0:nfft/2-1)*Fs/nfft;

loglog(f,mx);

Ideally I should be getting a single peak which corresponds to the resonance/natural frequency of the cantilever. Surprisingly, I am getting multiple peaks with the peaks located at different frequencies for different dataset (some peaks coincide). This is due to noise from other sources. How can I get rid of the undesired frequencies? Any algorithm to filter out the noises? Is there an alternate way to analyse the dataset to obtain the natural frequency?

Regards,

Rajesh

[Image Analyst]

How do you know there is noise from other sources? Why can't the cantilever have a variety of frequencies present? Are you moving it with a perfect sine wave? If not then there may be multiple spikes, or not even any spikes at all depending on how it moves.

[Rajesh]

Currently I am observing the fluctuation in the cantilever position with no external force in effect. In the absence of an externally driving force the cantilever should vibrate at its own natural frequency. I am expecting noise from other sources because of mechanical and structural vibrations (other parts of the set up, vibration from the building etc.) present while I am doing the experiment. The cantilever I am studying (micro capillary) is quite sensitive to noise.

By variety of frequencies do you mean different modes of resonance? Yes that is possible but the frequencies should be integer multiples of each other. I don't see such a correlation.

I am observing multiple peaks with amplitude of oscillations in the range of 100s of microns. How do I figure out what is relevant and what isn't?

[Image Analyst]

Hmmm... your middle paragraph tells me that you don't have a good concept of Fourier (spectrum) analysis, and don't really know what the Fourier Transform means or how to interpret it. That's fine though, not everyone had a full semester course in it, like I did in grad school. But it would be good for you to at least learn the basics. I suggest you Google around and look for Fourier tutorials, especially those that talk about periodic signals.

Answer 1

[Star Strider]

This is not a trivial problem. I took this approach and got the result in the plot. You will have to experiment to get the result you want, likely adding one or more additional sine terms to get a better fit. You will need the Signal Processing Toolbox to design and implement the filter. I used the core MATLAB function fminsearch. My objective function will work with the other curve fitting functions in the Statistics Toolbox and Optimization Toolbox. The fundamental natural frequency is 1/b(3) Hz. This should get you started, although you will need to experiment to get the result you want.

The code:

data=xlsread('Rajesh sample_data_2.xlsx');
t=data(:,1);
y=data(:,3);
Fs=500;
Fn = Fs/2;                                                                      % Nyquist Frequency
L = size(t,1);
y = detrend(y);                                                                 % Remove Linear Trend

fty = fft(y)/L;
Fv = linspace(0, 1, fix(L/2)+1)*Fn;                                             % Frequency Vector
Iv = 1:length(Fv);                                                              % Index Vector

figure(1)
plot(Fv, abs(fty(Iv)))
grid

Wp =   5/Fn;                                                                    % Passband (Normalised)
Ws =  20/Fn;                                                                    % Stopband (Normalised)
Rp =  1;                                                                        % Passband Ripple
Rs = 25;                                                                        % Stopband Ripple
[n,Wn]  = buttord(Wp,Ws,Rp,Rs);                                                 % Filter Order
[b,a]   = butter(n,Wn);                                                         % Transfer Function Coefficients
[sos,g] = tf2sos(b,a);                                                          % Second-Order-Section For Stability

figure(2)
freqz(sos, 2048, Fs)

yf = filtfilt(sos,g,y);

yu = max(yf);
yl = min(yf);
yr = (yu-yl);                                                                   % Range of ‘y’
yz = yf-yu+(yr/2);
zt = t(yz(:) .* circshift(yz(:),[1 0]) <= 0);                                   % Find zero-crossings
per = 2*mean(diff(zt));                                                         % Estimate period
ym = mean(y);                                                                   % Estimate offset

fit = @(b,x)  b(1) .* exp(b(2).*x) .* (sin(2*pi*x./b(3) + 2*pi/b(4))) + b(5);   % Objective Function to fit
fcn = @(b) sum((fit(b,t) - yf).^2);                                             % Least-Squares cost function
s = fminsearch(fcn, [yr; -1;  per;  -1;  ym])                                   % Minimise Least-Squares

xp = linspace(min(t),max(t));

figure(3)
plot(t,y,'b')
hold on
plot(t,yf,'y', 'LineWidth',2)
plot(xp,fit(s,xp), 'r', 'LineWidth',1.5)
hold off
grid
title('Sample Data 2')
xlabel('Time')
ylabel('Amplitude')
legend('Original Data', 'Lowpass-Filtered Data', 'Fitted Curve')

[uzzi]

Hello,

I was also trying this code. It is working perfectly well until the line

s = fminsearch(fcn, [yr; -1;  per;  -1;  ym])

where the error message is showing as

I also tried as

s = fminsearch(fcn, [double(yr); -1;  per;  -1;  double(ym)])

but the error is still there.

Can you help me with this, please, @Star Strider?

[Star Strider]

@uzzi — Please see my Comment.

Answer 2

[Image Analyst]

Depending on how it's moving, You may or may not see spikes. From what you say, I see no reason to expect it have spikes unless the structure is vibrating. If it's a sinusoid (in "resonance"), you'll see two spikes. If there are "harmonics" you will see multiple spikes as spacings a multiple of the first harmonic. But in general you will probably see a mountain shaped spectrum (due to the non-sinusoidal, wild and crazy motions that are present) with perhaps some small peaks on the sides of the mountain but only if it's vibrating. In other words, it is what it is. You got what you got. If you don't see an array of spikes, then that's just the way it is and you'll have to deal with it.