Advertisement
Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- a=np.random.normal(0,1,1024)
- f1=500
- f2=1200
- Fs = 8000
- x = np.arange(1024)
- b = 0.5*np.sin(2 * np.pi * f1 * x / Fs)+np.sin(2 * np.pi * f2 * x / Fs)
- c=b+0.1*a
- plt.figure(figsize=(15,4))
- fax, Pxx =signal.periodogram(a)
- Pxx[0]=0
- plt.plot(fax,10*np.log10(abs(Pxx)))
- plt.ylabel('Power/frequency (dB/rad/sample)')
- plt.xlabel('Normalized Frequency (xpirad/sample)')
- plt.title('Periodogram Power Spectral Densinity Estimate',fontweight='bold')
- r=[0.0,0.1,0.2,0.3,0.4,0.5]
- plt.xticks(r,('0','0.2','0.4','0.6','0.8','1'))
- plt.show()
- plt.figure(figsize=(15,4))
- fax, Pxx =signal.periodogram(b)
- Pxx[0]=0
- plt.plot(fax,10*np.log10(abs(Pxx)))
- plt.ylabel('Power/frequency (dB/rad/sample)')
- plt.xlabel('Normalized Frequency (xpirad/sample)')
- plt.title('Periodogram Power Spectral Densinity Estimate',fontweight='bold')
- r=[0.0,0.1,0.2,0.3,0.4,0.5]
- plt.xticks(r,('0','0.2','0.4','0.6','0.8','1'))
- plt.show()
- plt.figure(figsize=(15,4))
- fax, Pxx =signal.periodogram(c)
- Pxx[0]=0
- plt.plot(fax,10*np.log10(abs(Pxx)))
- plt.ylabel('Power/frequency (dB/rad/sample)')
- plt.xlabel('Normalized Frequency (xpirad/sample)')
- plt.title('Periodogram Power Spectral Densinity Estimate',fontweight='bold')
- r=[0.0,0.1,0.2,0.3,0.4,0.5]
- plt.xticks(r,('0','0.2','0.4','0.6','0.8','1'))
- plt.show()
- plt.figure(figsize=(15,4))
- fax, Pxx =signal.periodogram(signal.hann(1024)*a)
- Pxx[0]=0
- plt.plot(fax,10*np.log10(abs(Pxx)))
- plt.ylabel('Power/frequency (dB/rad/sample)')
- plt.xlabel('Normalized Frequency (xpirad/sample)')
- plt.title('Periodogram Power Spectral Densinity Estimate',fontweight='bold')
- r=[0.0,0.1,0.2,0.3,0.4,0.5]
- plt.xticks(r,('0','0.2','0.4','0.6','0.8','1'))
- plt.show()
- plt.figure(figsize=(15,4))
- fax, Pxx =signal.periodogram(signal.hann(1024)*b)
- Pxx[0]=0
- plt.plot(fax,10*np.log10(abs(Pxx)))
- plt.ylabel('Power/frequency (dB/rad/sample)')
- plt.xlabel('Normalized Frequency (xpirad/sample)')
- plt.title('Periodogram Power Spectral Densinity Estimate',fontweight='bold')
- r=[0.0,0.1,0.2,0.3,0.4,0.5]
- plt.xticks(r,('0','0.2','0.4','0.6','0.8','1'))
- plt.show()
- plt.figure(figsize=(15,4))
- fax, Pxx =signal.periodogram(signal.hann(1024)*c)
- Pxx[0]=0
- plt.plot(fax,10*np.log10(abs(Pxx)))
- plt.ylabel('Power/frequency (dB/rad/sample)')
- plt.xlabel('Normalized Frequency (xpirad/sample)')
- plt.title('Periodogram Power Spectral Densinity Estimate',fontweight='bold')
- r=[0.0,0.1,0.2,0.3,0.4,0.5]
- plt.xticks(r,('0','0.2','0.4','0.6','0.8','1'))
- plt.show()
- plt.figure(figsize=(9,7))
- fax, Pxx =signal.welch(a*signal.hann(1024),8000,nfft=256,nperseg=128)
- Pxx[0]=0
- plt.plot(fax,10*np.log10(abs(Pxx)),label='NFFT=256')
- fax, Pxx =signal.welch(a*signal.hann(1024),8000,nfft=128,nperseg=64)
- Pxx[0]=0
- plt.plot(fax,10*np.log10(abs(Pxx)),label='NFFT=128')
- fax, Pxx =signal.welch(a*signal.hann(1024),8000,nfft=64,nperseg=32)
- Pxx[0]=0
- plt.plot(fax,10*np.log10(abs(Pxx)),label='NFFT=64')
- plt.legend()
- plt.ylabel('PSD (dB/Hz)')
- plt.xlabel('Frequency (Hz)')
- plt.title('Szum gaussowski - Periodogram Welscha',fontweight='bold')
- plt.figure(figsize=(9,7))
- fax, Pxx =signal.welch(b*signal.hann(1024),8000,nfft=256,nperseg=128)
- Pxx[0]=0
- plt.plot(fax,10*np.log10(abs(Pxx)),label='NFFT=256')
- fax, Pxx =signal.welch(b*signal.hann(1024),8000,nfft=128,nperseg=64)
- Pxx[0]=0
- plt.plot(fax,10*np.log10(abs(Pxx)),label='NFFT=128')
- fax, Pxx =signal.welch(b*signal.hann(1024),8000,nfft=64,nperseg=32)
- Pxx[0]=0
- plt.plot(fax,10*np.log10(abs(Pxx)),label='NFFT=64')
- plt.legend()
- plt.ylabel('PSD (dB/Hz)')
- plt.xlabel('Frequency (Hz)')
- plt.title('Szum gaussowski - Periodogram Welscha',fontweight='bold')
- plt.figure(figsize=(9,7))
- fax, Pxx =signal.welch(c*signal.hann(1024),8000,nfft=256,nperseg=128)
- Pxx[0]=0
- plt.plot(fax,10*np.log10(abs(Pxx)),label='NFFT=256')
- fax, Pxx =signal.welch(c*signal.hann(1024),8000,nfft=128,nperseg=64)
- Pxx[0]=0
- plt.plot(fax,10*np.log10(abs(Pxx)),label='NFFT=128')
- fax, Pxx =signal.welch(c*signal.hann(1024),8000,nfft=64,nperseg=32)
- Pxx[0]=0
- plt.plot(fax,10*np.log10(abs(Pxx)),label='NFFT=64')
- plt.legend()
- plt.ylabel('PSD (dB/Hz)')
- plt.xlabel('Frequency (Hz)')
- plt.title('Szum gaussowski - Periodogram Welscha',fontweight='bold')
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement