一个简单的python数字锁相环

一个简单的数字锁相环demo：

import matplotlibmatplotlib.use('Qt5Agg')import numpy as npimport matplotlib.pyplot as plt#parametersphase_offset      = 0.00    #carrier phase offsetfrequency_offset  = 0.30    #carrier frequency offsetwn                = 0.01    #pll bandwidthzeta              = 0.707   #pll damping factorK                 = 1000    # pll loop gainn                 = 400   #number of samples#generate loop filter parameters (active PI design)t1 = K/(wn*wn)    #tau_1t2 = 2*zeta/wn    #tau_2#feed-forward coefficients (numerator)b0 = (4*K/t1)*(1.+t2/2.0)b1 = (8*K/t1)b2 = (4*K/t1)*(1.-t2/2.0)#feed-back coefficients (denominator)a1 = -2.0a2 =  1.0#print filter coefficients (as comments)print b0,b1,b2print a1,a2#filter bufferv0=0.0v1=0.0v2=0.0#initialize statesphi     = phase_offset   #input signal's initial phasephi_hat = 0.0            #PLL's initial phasen =500delta_phi = np.zeros(n)#run basic simulationfor i in range(n):    #compute input sinusoid and update phase    x = np.exp(1j*phi)    phi += frequency_offset    #compute PLL output from phase estimate    y = np.exp(1j * phi_hat)    #compute error estimate    delta_phi[i] = np.angle(x * np.conjugate(y))#cargf( x * conjf(y) )carg is a standard library function which calculates the argument (phase angle) of a complex number.    #push result through loop filter, updating phase estimate    #advance buffer    v2 = v1  #shift center register to upper register    v1 = v0  #shift lower register to center register    #compute new lower register    v0 = delta_phi[i] - v1*a1 - v2*a2    #compute new output    phi_hat = v0*b0 + v1*b1 + v2*b2    phi_hat = np.mod(phi_hat, 2 * np.pi)print delta_phiplt.plot(delta_phi)plt.grid()plt.show()

可以看出迭代近200次后收敛。