I need to find power spectrum of a signal applying Welch method. Then using Chebyshev, Kaiser and Gauss windows, find the best window in terms of frequency separation.
I've wrote a function that produces noised sinusoidal signal and used the Welch method to this signal, also generated all of these three window functions.
import numpy as np
import matplotlib.pyplot as plt
import scipy.signal as sg
def sig_noise(f,snr,n):
fs=1000
dt=1/fs
tk=n*dt
A=1
B=10**(-snr/20)
t=np.arange(0,tk,dt)
y=(np.random.rand(n)-0.5)
z=A*np.sin(2*np.pi*f*t)+B*y
return z,t
[s,t]=sig_noise(140,-10,64)
f, Pxx_den = sg.welch(s,nperseg=64,noverlap=None) # Welch method, parameters are taken from book
plt.plot(f,Pxx_den)
plt.xlabel("Frequency")
plt.grid()
plt.figure(figsize=(10,10))
chebyshev=sg.chebwin=(64,40) # Chebyshev window
kaiser=sg.kaiser(64,0) # Kaiser window
gauss=sg.gaussian(64,np.std((f,Pxx_den),ddof=0)) # Gauss window
# plots for check
plt.subplot(221)
plt.grid()
plt.title("Gauss")
plt.plot(gauss)
plt.subplot(222)
plt.grid()
plt.title("Chebyshev")
plt.plot(chebyshev)
plt.subplot(223)
plt.grid()
plt.title("Kaiser")
plt.plot(kaiser)
f, Pxx_den = sg.welch(s,nperseg=64,noverlap=None,window=chebyshev)
the last line of code doesn't work as I want to use for example chebyshev window. Unknown window type is what I get as the answer. Could somebody tell me where is the mistake?
The issue is with a extra =
character on the following line:
chebyshev=sg.chebwin=(64,40) # Chebyshev window
# ^^^
# extra "=" character
This causes chebyshev
to be equal to the tuple (64,40)
rather than the intended array (following a similar approach as for your other windows). Later when given a tuple as windows argument, scipy.signal.welch
tries to construct the window by calling scipy.signal.get_window
for a window type specified by the first element of the tuple. Obviously 64
(the first element of the tuple in your case) is not a valid window type.
So, the fix is to simply correct the typo:
chebyshev=sg.chebwin(64,40) # Chebyshev window
or alternatively call welch
with an appropriately constructed tuple:
f, Pxx_den = sg.welch(s,nperseg=64,noverlap=None,window=('chebwin',40)