Finding the dominant tone in a signal

I am given 256 samples at a sampling rate of 48k and am asked to find the frequency of dominate tones and the amplitude of the largest tone. The data looks like this when plotted:

I assume I must use a FFT transform to do this. I've tried this in matlab:

Fs = 48000;
x = dataSet; % 256 samples data set
xdft = fft(x);
maxAmp = max(abs(xdft));
%Not sure how I can grab the tones here
freq = 0:Fs / length(x):Fs/2;


Not really sure how to proceed, Any ideas?

I think you are having trouble constructing your frequency axis properly. Once you do this, you can do a simple peak pick. I have re-written the code for you:

Fs = 48000;
x = dataSet;            % 256 samples data set
fftLength = length(x);  % Always make sure to be at least as long as your data.
xdft = fft(x,fftLength);
maxAmp = max(abs(xdft));

freq = [0:fftLength-1].*(Fs/fftLength); % This is your total freq-axis
freqsYouCareAbout = freq(freq < Fs/2);  % You only care about either the pos or neg
% frequencies, since they are redundant for
% a real signal.

xdftYouCareAbout = abs(xdft(1:round(fftLength/2))); % Take the absolute magnitude.

maxFreq = freqsYouCareAbout(index); % This is the frequency of your dominant signal.


Have you tried plotting the magnitude of xdft? You should be able to see where the peaks are. To find the frequency associated with any peak, find the index where the peak magnitude lives, and use that to pull the frequency from an array like your freq, except you should create it from 0 to Fs with

freq = 0:Fs/(length(xdft)-1):Fs;


since the DFT runs from 0 to the sampling frequency. If you would prefer it in the range -Fs/2:Fs/2, you can reorder it with fftshift, but I'll let you read up on that function. You should also find that the spectrum is symmetrical, so you'll get two peaks (of the same height) for the strongest frequency. You can easily drop half of the spectrum for processing, of course.

OK, there are several issues here:

1. You need to detrend :-) the data to remove the non-zero mean.

x = dataSet - mean(dataSet);

This is because the largest component (when I do an FFT) is the DC offset.

2. Then you can use any number of techniques to find the next highest. The previous answers are good. This link has an implementation of it that I did many moons ago called discperiod.

omega1 = discperiod(dataSet);

3. That, however, will only find the next highest peak in the spectrum. So you need to remove the frequency component found in step 2. To do that in a statistically sensible way, you need to find the amplitude and phase of that component. One way to do this is to solve the least squares problem: $$\min_{A, \phi} \sum_{t=0}^{T-1} \left | x[t] - A \cos({\tt omega1} t + \phi) \right|^2$$ That is, minimize over $A$ and $\phi$. This can be done pretty simply by just finding the FFT coefficient of omega1. Julius O. Smith III derives this here. The mean correction is then

x_res = x - A*cos(omega*[0:T-1]'+phi);

4. Trying to find the next highest peak from x_res can show up a problem: the next highest peak is probably on the "edge" of the first peak (that has just been removed). This is because the discperiod function quantizes the possible frequencies to the FFT bin frequencies. A simple way to fix this is to choose what Mohammad calls fftLength to be as large as the accuracy required. This can be several times larger than the length of the signal.

No more time! Will try to write up more later.