from pyestimate.estimators import sin_param_estimate
import numpy as np
import matplotlib.pyplot as plt
N = 20 # number of samples
f = 0.123456 # frequency to be estimated
A = 1.23456 # amplitude to be estimated
phi = np.pi/7 # phase to be estimated
sigma = 1 # standard deviation of WGN
n = np.arange(N)
s = A * np.cos(2*np.pi*f*n+phi) # original signal
w = np.random.default_rng(seed=0).normal(scale=sigma, size=N) # white gaussian noise
x = s+w # input signal for estimation: sine + noise
#
# Estimate sinusoid parameters
#
A_hat, f_hat, phi_hat = sin_param_estimate(x) # parameters estimation
#
# Reconstruct original signal from estimated parameters
#
s_hat = A_hat * np.cos(2*np.pi*f_hat*n+phi_hat) # estimated signal
#
# Plot the original signal, the input signal corrupted with noise and the reconstructed signal
#
plt.plot(n, s, linewidth=3.0, label='original signal')
plt.plot(n, x, label='corrupted signal')
plt.plot(n, s_hat, 'k--', label='estimated signal')
plt.xlabel('$n$')
plt.ylabel('$x[n]$')
plt.title('Sinusoidal frequency, amplitude and phase estimation in WGN')
plt.legend()
plt.grid()
plt.show()