FFT and Power spectral estimation
- Power Spectrum---Periodic graph method for signal extraction using Fourier transform
ClF
fs=1000;
n=256; The length of the nfft=256;% data and the length of the data used by the FFT
Time series used by the n=0:n-1;t=n/fs;%
Xn=sin (2*pi*50*t) +2*sin (2*pi*120*t) +randn (1,n);
PXX=10*LOG10 (ABS (FFT (XN,NFFT). ^2)/N),%fourier average of the amplitude spectrum squared and converted to DB
F= (0:length (PXX)-1) *fs/length (Pxx);% gives frequency sequence
Subplot (2,1,1), Plot (F,PXX),% plotting power spectral curves
Xlabel (' frequency/hz '); Ylabel (' Power spectrum/db ');
Title (' Cycle graph n=256 '); grid on;
fs=1000;
n=1024; The length of the nfft=1024;% data and the length of the data used by the FFT
Time series used by the n=0:n-1;t=n/fs;%
Xn=sin (2*pi*50*t) +2*sin (2*pi*120*t) +randn (1,n);
PXX=10*LOG10 (ABS (FFT (XN,NFFT). ^2)/N),%fourier average of the amplitude spectrum squared and converted to DB
F= (0:length (PXX)-1) *fs/length (Pxx);% gives frequency sequence
Subplot (2,1,2), Plot (F,PXX),% plotting power spectral curves
Xlabel (' frequency/hz '); Ylabel (' Power spectrum/db ');
Title (' Cycle graph n=256 '); grid on;
- Power Spectrum---Segmented periodic graph method for signal extraction using Fourier transform
% thought: dividing the signal into overlapping or non-overlapping segments, power spectral estimation of each small sequence of signals, then averaging as the power spectrum of the entire sequence
ClF
fs=1000;
n=1024; The length of the nsec=256;% data and the length of the data used by the FFT
Time series used by the n=0:n-1;t=n/fs;%
RANDN (' state ', 0);
Xn=sin (2*pi*50*t) +2*sin (2*pi*120*t) +randn (1,n);
Pxx1=abs (FFT (xn (1:256), Nsec). ^2)/nsec; % first segment Power spectrum
Pxx2=abs (FFT (xn (257:512), Nsec) ^2)/nsec;% second segment power spectrum
Pxx3=abs (FFT (xn (513:768), Nsec). ^2)/nsec;% third segment power spectrum
Pxx4=abs (FFT (xn (769:1024), Nsec) ^2)/nsec;% fourth segment power spectrum
PXX=10*LOG10 (PXX1+PXX2+PXX3+PXX4/4);%fourier average of amplitude spectrum squared and converted to DB
F= (0:length (PXX)-1) *fs/length (Pxx);% gives frequency sequence
Subplot (2,1,1), Plot (f (1:NSEC/2), Pxx (1:NSEC/2)),% plot power spectral curve
Xlabel (' frequency/hz '); Ylabel (' Power spectrum/db ');
Title (' Average period graph (no overlap) n=4*256 '); grid on;
% estimation of power spectrum using signal overlap segment
Pxx1=abs (FFT (xn (1:256), Nsec). ^2)/nsec; % first segment Power spectrum
Pxx2=abs (FFT (xn (129:384), Nsec) ^2)/nsec;% second segment power spectrum
Pxx3=abs (FFT (xn (257:512), Nsec). ^2)/nsec;% third segment power spectrum
Pxx4=abs (FFT (xn (385:640), Nsec) ^2)/nsec;% fourth segment power spectrum
Pxx5=abs (FFT (xn (513:768), Nsec) ^2)/nsec;% fourth segment power spectrum
Pxx6=abs (FFT (xn (641:896), Nsec) ^2)/nsec;% fourth segment power spectrum
Pxx7=abs (FFT (xn (769:1024), Nsec) ^2)/nsec;% fourth segment power spectrum
PXX=10*LOG10 (PXX1+PXX2+PXX3+PXX4+PXX5+PXX6+PXX7/7);%fourier average of amplitude spectrum squared and converted to DB
F= (0:length (PXX)-1) *fs/length (Pxx);% gives frequency sequence
Subplot (2,1,2), Plot (f (1:NSEC/2), Pxx (1:NSEC/2)),% plot power spectral curve
Xlabel (' frequency/hz '); Ylabel (' Power spectrum/db ');
Title (' Average period graph (overlapping) n=1024 '); grid on;
- A power Spectrum---Welch method for finding signals using Fourier transform
% thought: The Welch method uses the signal overlapping segment, the window function and the FFT algorithm to calculate the self power spectrum of a signal sequence (PSD) and two signal sequence of the cross Power Spectrum (CSD), using MATLAB self-
% function PSD with
ClF
fs=1000;
n=1024; nfft=256;n=0:n-1;t=n/fs;
Window=hanning (256);
noverlap=128;
Dflag= ' None ';
RANDN (' state ', 0);
Xn=sin (2*pi*50*t) +2*sin (2*pi*120*t) +randn (1,n);
PXX=PSD (Xn,nfft,fs,window,noverlap,dflag);
f= (0:NFFT/2) *fs/nfft;
Plot (F,10*LOG10 (Pxx));
Xlabel (' frequency/hz '); Ylabel (' Power spectrum/db ');
Title (' Psd--welch method '); grid on;
- Power Spectrum estimation----Multi-Window Method (Multitaper method, MTM)
% idea: Using multiple orthogonal windows to obtain their own independent approximate power spectral estimation, synthesize these to obtain a sequence of power spectral estimation; more freedom than normal periodic graphs The MTM method takes one parameter: Time with a
% wide product NW, which defines the number of Windows used to calculate the power spectrum as 2*NW-1,NW, the higher the time domain resolution and the lower the frequency resolution, the less the fluctuation of the power spectral estimate; with the NW increase
%, the spectral leakage in each estimate increases , the deviation of the total power spectral estimation is increased by
CLF;
fs=1000;
n=1024; nfft=256;n=0:n-1;t=n/fs;
Randn (' state ', 0);
Xn=sin (2*pi*50*t) +2*sin (2*pi*120*t) +randn (1,n);
[Pxx1,f]=pmtm (XN,4,NFFT,FS);% There is a problem
subplot (2,1,1), Plot (F,10*log10 (Pxx1));
Xlabel (' frequency/hz '); Ylabel (' Power spectrum/db ');
Title (' Multi-Window method (MTM) nw=4 '); grid on;
[Pxx,f]=pmtm (XN,2,NFFT,FS);
Subplot (2,1,2), Plot (F,10*log10 (Pxx));
Xlabel (' Frequency/hz '); Ylabel (' Power spectrum/db ');
Title (' Multi-Window method (MTM) nw=2 '); grid on;
- Power spectrum Estimation----Maximum entropy method (Maxmum entmpy Method,mem method)
% thought: the assumption that the random sequence is a stationary Gaussian process utilizes known autocorrelation sequences rxx (0), Rxx (1), Rxx (2) ... rxx (p) as the basis, Extrapolation autocorrelation sequence Rxx (p+1), Rxx (p+2) ... Guaranteed maximum information entropy
CLF;
fs=1000;
n=1024; nfft=256;n=0:n-1;t=n/fs;
Window=hanning (256);
Randn (' state ', 0);
Xn=sin (2*pi*50*t) +2*sin (2*pi*120*t) +randn (1,n);
[Pxx1,f]=pmem (XN,14,NFFT,FS);% There is a problem
subplot (2,1,1), Plot (F,10*log10 (Pxx1));
Xlabel (' frequency/hz '); Ylabel (' Power spectrum/db ');
Title (' Maximum entropy method (MEM) order=14 '); Gridon; The
% uses the Welch method to estimate the power spectrum
noverlap=128;
Dflag= ' None ';
Subplot (2,1,2)
PSD (xn,nfft,fs,window,noverlap,dflag);
Xlabel (' frequency/hz '); Ylabel (' Power spectrum/db ');
Title (' Welch Method estimates power spectrum '); grid on;
- Power Spectrum estimation----Multi-signal classification (multiple signal Classification,music method)
% Note: For multi-sine frequency estimation in white noise
% thought: The data autocorrelation matrix is composed of the signal autocorrelation matrix and the noise autocorrelation matrix, and the matrix eigenvalue vector is obtained.
ClF
fs=1000;
n=1024; nfft=256;n=0:n-1;t=n/fs;
RANDN (' state ', 0);
Xn=sin (2*pi*100*t) +2*sin (2*pi*200*t) +randn (1,n);
Pmusic (xn,[7,1.1],nfft,fs,32,16);
Xlabel (' frequency/khz '); Ylabel (' Power spectrum/db ');
Title (' Welch Method estimates power spectrum '); grid on;
Matlab for FFT and Power spectrum
