1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24
|
import matplotlib.pyplot as plt
import numpy as np
import scipy.signal as signal
rng = np.random.default_rng(73625) # seeding for reproducibility
fs, n = 10e3, 10_000
f_x, noise_power = 1270, 1e-3 * fs / 2
t = np.arange(n) / fs
x = (np.sqrt(2) * np.sin(2 * np.pi * f_x * t) +
rng.normal(scale=np.sqrt(noise_power), size=t.shape))
fg, axx = plt.subplots(1, 2, sharex='all', tight_layout=True, figsize=(7, 3.5))
axx[0].set(title="Squared Magnitude Spectrum", ylabel="Square of Magnitude in V²")
axx[1].set(title="Power Spectral Density", ylabel="Power Spectral Density in V²/Hz")
for ax_, s_ in zip(axx, ('spectrum', 'density')):
f_p, P_p = signal.periodogram(x, fs, 'hann', scaling=s_)
f_w, P_w = signal.welch(x, fs, scaling=s_)
ax_.semilogy(f_p/1e3, P_p, label=f"Periodogram ({len(f_p)} bins)")
ax_.semilogy(f_w/1e3, P_w, label=f"Welch's Method ({len(f_w)} bins)")
ax_.set(xlabel="Frequency in kHz", xlim=(0, 2), ylim=(1e-7, 1.3))
ax_.grid(True)
ax_.legend(loc='lower center')
plt.show()
|