def simple_SA(x,NS,NFFT,fs,NAVG=<NUM_LIT:1>,window='<STR_LIT>'):
Nx = len(x)<EOL>K = int(Nx/NS)<EOL>print('<STR_LIT>', K)<EOL>if NAVG > K:<EOL><INDENT>print('<STR_LIT>')<EOL>return <NUM_LIT:0>,<NUM_LIT:0><EOL><DEDENT>if window.lower() == '<STR_LIT>' or window.lower() == '<STR_LIT>':<EOL><INDENT>w = signal.boxcar(NS)<EOL><DEDENT>elif window.lower() == '<STR_LIT>':<EOL><INDENT>w = signal.hanning(NS)<EOL><DEDENT>xsw = np.zeros((K,NS)) + <NUM_LIT>*np.zeros((K,NS))<EOL>for k in range(NAVG):<EOL><INDENT>xsw[k,] = w*x[k*NS:(k+<NUM_LIT:1>)*NS]<EOL><DEDENT>Sx = np.zeros(NFFT)<EOL>for k in range(NAVG):<EOL><INDENT>X = fft.fft(xsw[k,],NFFT)<EOL>Sx += abs(X)**<NUM_LIT:2><EOL><DEDENT>Sx /= float(NAVG)<EOL>Sx /= float(NFFT**<NUM_LIT:2>)<EOL>NFFTby2 = int(NFFT/<NUM_LIT:2>)<EOL>if x.dtype != '<STR_LIT>':<EOL><INDENT>n = np.arange(NFFTby2)<EOL>f = fs*n/float(NFFT)<EOL>Sx = Sx[<NUM_LIT:0>:NFFTby2]<EOL><DEDENT>else:<EOL><INDENT>n = np.arange(NFFTby2)<EOL>f = fs*np.hstack((np.arange(-NFFTby2,<NUM_LIT:0>),np.arange(NFFTby2)))/float(NFFT)<EOL>Sx = np.hstack((Sx[NFFTby2:],Sx[<NUM_LIT:0>:NFFTby2]))<EOL><DEDENT>return f, Sx<EOL>
Spectral estimation using windowing and averaging.

This function implements averaged periodogram spectral estimation
estimation similar to the NumPy's psd() function, but more
specialized for the the windowing case study of Chapter 16.

Parameters
----------
x : ndarray containing the input signal
NS : The subrecord length less zero padding, e.g. NS < NFFT
NFFT : FFT length, e.g., 1024 = 2**10
fs : sampling rate in Hz
NAVG : the number of averages, e.g., 1 for deterministic signals
window : hardcoded window 'boxcar' (default) or 'hanning'

Returns
-------
f : ndarray frequency axis in Hz on [0, fs/2]
Sx : ndarray the power spectrum estimate

Notes
-----
The function also prints the maximum number of averages K possible
for the input data record.

Examples
--------
>>> import matplotlib.pyplot as plt
>>> import numpy as np
>>> from sk_dsp_comm import sigsys as ss
>>> n = np.arange(0,2048)
>>> x = np.cos(2*np.pi*1000/10000*n) + 0.01*np.cos(2*np.pi*3000/10000*n)
>>> f, Sx = ss.simple_SA(x,128,512,10000)
>>> plt.plot(f, 10*np.log10(Sx))
>>> plt.ylim([-80, 0])
>>> plt.xlabel("Frequency (Hz)")
>>> plt.ylabel("Power Spectral Density (dB)")
>>> plt.show()

With a hanning window.

>>> f, Sx = ss.simple_SA(x,256,1024,10000,window='hanning')
>>> plt.plot(f, 10*np.log10(Sx))
>>> plt.xlabel("Frequency (Hz)")
>>> plt.ylabel("Power Spectral Density (dB)")
>>> plt.ylim([-80, 0])