MATLAB Examples

# Example 8.6.5. A Comparison of Frequency Estimation Methods.

## Workspace Initialization.

```clc; clear; close all; ```

## Signal Definition

```N = [64 128 256 512]; % Number of samples used for each case. NumOfRuns = 10; % Number of frequency estimation functions to be calculated for every different data record. % Open a new figure: h1 = figure('NumberTitle', 'off','Name', ... 'Figure 8.37 Frequency estimation for a process consisting of four complex exponentials in WGN', ... 'Visible','off','Position', [100 0 800 950]); sigmaw = sqrt(0.5); % Create frequency axis: L = 1024; w = 0:2*pi/L:2*pi*(1-1/L); Pphd = zeros(NumOfRuns,L); for m=1:length(N) for k=1:NumOfRuns ```
``` noise = sigmaw*randn(1,N(m)); n = 0:N(m)-1; phi = 2*pi*(rand(4,1) - 0.5); omegas = [0.2*pi 0.3*pi 0.8*pi 1.2*pi].'; % Compute N(m) samples of the given Harmonic Process: x = sum(exp(1i*omegas*n + phi*ones(1,N(m)))) + noise; ```

## Estimate the Power Spectrum using the Pisarenko Harmonic Decomposition (PHD):

``` p = 4; % Number of complex exponentials present assumed known. [vmin, sigma] = phd(x,p); for k1=1:L e = exp(-1i*w(k1)*(0:p)); Pphd(k,k1) = 1./abs(e*vmin).^2; end Pphd(k,:) = 10*log10(Pphd(k,:)); ```
``` end subplot(4,1,m) plot(w/pi,Pphd,'k','LineWidth',0.1) title(['Overlay of ',num2str(NumOfRuns),' PHD Spectra Using',num2str(N(m)),' samples']) grid on; axis tight; ylim([-10 100]) xlabel('Frequency (units of pi)'); end % Restore the visibility of the figure. set(h1, 'Visible','on'); ```