MATLAB Examples

Estimate the Power Spectrum Using dsp.SpectrumEstimator

Alternately, you can compute the power spectrum of the signal using the dsp.SpectrumEstimator System object. You can acquire the output of the spectrum estimator and store the data for further processing. To view other objects in the Estimation library, type help dsp in the MATLAB® command prompt, and click Estimation.

Contents

Initialization

Use the same source as in the previous section on using the dsp.SpectrumAnalyzer to estimate the power spectrum. The input sine wave has two frequencies: one at 1000 Hz and the other at 5000 Hz. Initialize dsp.SpectrumEstimator to compute the power spectrum of the signal using the filter bank approach. View the power spectrum of the signal using the dsp.ArrayPlot object.

Fs = 44100;
Sineobject1 = dsp.SineWave('SamplesPerFrame',1024,'PhaseOffset',10,...
    'SampleRate',Fs,'Frequency',1000);
Sineobject2 = dsp.SineWave('SamplesPerFrame',1024,...
    'SampleRate',Fs,'Frequency',5000);

SpecEst = dsp.SpectrumEstimator('Method','Filter bank',...
    'PowerUnits','dBm','SampleRate',Fs,'FrequencyRange','onesided');
ArrPlot = dsp.ArrayPlot('PlotType','Line','ChannelNames',{'Power spectrum of the input'},...
    'YLimits',[-80 30],'XLabel','Number of samples per frame','YLabel',...
    'Power (dBm)','Title','One-sided power spectrum with respect to samples');

Estimation

Stream in and estimate the power spectrum of the signal. Construct a for-loop to run for 5000 iterations. In each iteration, stream in 1024 samples (one frame) of each sine wave and compute the power spectrum of each frame. Add Gaussian noise with mean at 0 and a standard deviation of 0.001 to the input signal.

for Iter = 1:5000
    Sinewave1 = Sineobject1();
    Sinewave2 = Sineobject2();
    Input = Sinewave1 + Sinewave2;
    NoisyInput = Input + 0.001*randn(1024,1);
    PSoutput = SpecEst(NoisyInput);
    ArrPlot(PSoutput);
end

Using the filter bank approach, the spectral estimate has a high resolution and the peaks are precise with no spectral leakage.

Convert x-axis to Represent Frequency

By default, the array plot shows the power spectral data with respect to the number of samples per frame. The number of points on the x-axis equals the length of the input frame. The spectrum analyzer plots the power spectral data with respect to frequency. For a one-sided spectrum, the frequency varies in the range [0 Fs/2]. For a two-sided spectrum, the frequency varies in the range [-Fs/2 Fs/2]. To convert the x-axis of the array plot from sample-based to frequency-based, do the following:

  • Click on the Configuration Properties icon.
  • For a one-sided spectrum - On Main tab, set Sample increment to $Fs/FrameLength$ and X-offset to 0.
  • For a two-sided spectrum - On Main tab, set Sample increment to $Fs/FrameLength$ and X-offset to $-Fs/2$.

In this example, the spectrum is one-sided and hence, the Sample increment and X-offset are set to 44100/1024 and 0, respectively. To specify the frequency in kHz, set the Sample increment to 44.1/1024.

ArrPlot.SampleIncrement = (Fs/1000)/1024;
ArrPlot.XLabel = 'Frequency (kHz)';
ArrPlot.Title = 'One-sided power spectrum with respect to frequency';

for Iter = 1:5000
    Sinewave1 = Sineobject1();
    Sinewave2 = Sineobject2();
    Input = Sinewave1 + Sinewave2;
    NoisyInput = Input + 0.001*randn(1024,1);
    PSoutput = SpecEst(NoisyInput);
    ArrPlot(PSoutput);
end

Live Processing

The output of the dsp.SpectrumEstimator object contains the spectral data and is available for further processing. The data can be processed in real-time or it can be stored in the workspace.