The goal of frequency estimation is, given N successive observations of a pure sinusoidal or exponential tone, contaminated by AWGN noise, to estimate as accurately as possible the frequency of the signal. Note that the two problems (real and complex cases) are not equivalent.
These scripts (see links at the right of this page) contain the SCILAB implementation of various estimators, based both on time series analysis and on DFT interpolation. You can also find in the attached PDF note (link at the right) a short description of the available estimators and the results I obtained by simulating them.
First, if not already done, install a recent version of Scilab on your computer (official SCILAB download page). Then, download the two scripts (links at the right of this page) and execute them inside Scilab (to do this, just open them with the Scilab editor Scinotes, and press F5). Now you can use the estimation functions as described below.
estimator = get_freqestim(name)
Build an estimator object. Supported estimators are the following (see attached PDF file for a short description of each one):
'kc', 'lrp', 'pscfd', 'argmax(dft)', 'quadratic',
'Gaussian', 'McEachern', 'Jacobsen', 'Candan',
'Candan 2', 'Cedron', 'Cedron cplx'.
Note that the name specification is not case sensitive.
nf = freqestim(estimator, x)
Estimate the normalized frequency nf (nf = f/fs), using the specified estimator and the input signal x.
Test a list of frequency estimators: plots the mean square error versus SNR for each estimator. lst must be a list of estimators objects, cols is a vector of color specification (same as plot third argument), and form specify real ('r') or complex ('c') test.
fs = 1e3; // 1 KHz sampling frequency f = 0.21 * fs; // F = 0.21 Fs = 210 Hz N = 128; // 128 samples // Get the Jacobsen estimator, designed for complex tones estimator = get_freqestim('Jacobsen'); // Generate the test signal t = linspace(0,(N-1)/fs,N); x = exp(2*%pi*%i*f*t); // Proceed to estimation f = freqestim(estimator, x) disp(f) 0.2100002
fs = 1e3; // 1 KHz sampling frequency f = 0.21 * fs; // F = 0.21 Fs = 210 Hz N = 128; // 128 samples // Get the Cedron estimator, designed for real tones estimator = get_freqestim('Cedron'); // Generate the test signal (real signal this time) t = linspace(0,(N-1)/fs,N); x = sin(2*%pi*f*t); // Proceed to estimation f = freqestim(estimator, x) disp(f) 0.21
This example shows how to use the function
test_freqestim so as to compare different algorithms (here Gaussian, Candan 2, Cedron and PSCFD estimators, for complex tones).
// Build the list of frequency estimators lst= list(get_freqestim('Gaussian'),get_freqestim('Candan 2'),... get_freqestim('Cedron'), get_freqestim('PSCFD')); // Curve colors: red+stars, cyan+diamonds, magenta, blue+circles cols = ['r-*','c-d','m-','b-o']; // Test and plot MSE = f(SNR), complex case test_freqestim(lst, cols, 'c')
Here is the figure obtained:
Below are more complete graphs, both for complex and real tones detection: