I have data from a sensor and I need to find the frequency of it. It looks like fft() seems to be the way to go, but the MATLAB docs only show how to get a graph of the frequencies, I don't know what to do from there.
Here's what my data looks like:
One way to go is indeed to use an fft. Since the fft gives you the frequency representation of the signal, you want to look for the maximum, and since the fft is a complex signal, you will want to take the absolute value first. The index will correspond to the normalized frequency with maximum energy. Last, if your signal has an offset, as is the case with the one you show, you want to get rid of that offset before taking the fft so that you do not get a max at the origin representing the DC component.
Everything I described put in one line would be:
[maxValue,indexMax] = max(abs(fft(signal-mean(signal))));
where indexMax is the index where the max fft value can be found.
Note: to get from indexMax to the actual frequency of interest, you will need to know the length L of the fft (same as the length of your signal), and the sampling frequency Fs. The signal frequency will then be:
frequency = indexMax * Fs / L;
Alternatively, faster and working fairly well too depending on the signal you have, take the autocorrelation of your signal:
autocorrelation = xcorr(signal);
and find the first maximum occurring after the center point of the autocorrelation. (The autocorrelation will be symmetric with its maximum in the middle.) By finding that maximum, you find the first place where the shifted signal looks more or less like itself. I.e. you find the period of your signal. Since the signal shifted by a multiple of its period will always look like itself, you need to make sure that the maximum you find indeed corresponds to the period of the signal and not one of its multiples.
Because of the noise in your signal, the absolute maximum could very well occur at a multiple of your period instead of the period itself. So to account for that noise, you would take the absolute max of the autocorrelation (autocorrelation(length(autocorrelation)/2+1), and then find where the autocorrelation is larger than, say, 95% of that maximum value for the first time in the second half of the signal. 95%, 99%, or some other number would depend on how much noise corrupts your signal.
UPDATE: I realize that I assumed you meant by "frequency" of your signal the pitch or base harmonic or frequency with the most energy, however you want to look at it. If by frequency you meant the frequency representation of your signal, then to a first approximation, you just want to plot the abs of the FFT to get an idea of where the energy is:
plot(abs(fft));
If you want to understand why there is an abs, or what relevant info you are losing by not representing the phase of the fft, you may want to read a bit more about the DFT transform to understand exactly what you get.
I think it should be
(indexMax-1) * Fs / L
The first element of abs(fft(x)) is the direct current (DC), or bias, or mean of signal, or X0. We count from the second element (X1).
Please let me know if I am wrong. Thanks.
clear all
clc
close all
Fs = 1;
T = 11 % Note this T is deliberately chosen , so that we have about 1.7 cycle of cosine singal
t = 0:Fs:T; % T seconds
L = length(t); % L is the length of sample sequence
bias = 4
signal = sin(t) + bias;
[maxValue,indexMax] = max(abs(fft(signal-mean(signal))));
frequency_method1 = (indexMax-1) * Fs / (L-1);
frequency_method2 = (indexMax-1) * Fs / L;
number_of_cycles_method1 = frequency_method1*T
number_of_cycles_method2 = frequency_method2*T
subplot(2,1,1)
plot(t,signal,'-or') ; grid on;
legend('about 1.7 cycles of cosine signal')
subplot(2,1,2)
plot(abs(fft(signal-mean(signal))),'-xb'); grid on
legend('abs of fft')
number_of_cycles_method1 =
2
number_of_cycles_method2 =
1.8333
