# NLMS算法的模擬測試 相關文檔： [_自適應濾波器和NLMS模擬_](fast_nlms_in_python.html) 測試NLMS在系統辨識、信號預測和信號均衡方面的應用。 ``` # -*- coding: utf-8 -*- # filename: nlms_test.py import numpy as np import pylab as pl import nlms_numpy import scipy.signal # 隨機產生FIR濾波器的系數，長度為length， 延時為delay， 指數衰減 def make_path(delay, length): path_length = length - delay h = np.zeros(length, np.float64) h[delay:] = np.random.standard_normal(path_length) * np.exp( np.linspace(0, -4, path_length) ) h /= np.sqrt(np.sum(h*h)) return h def plot_converge(y, u, label=""): size = len(u) avg_number = 200 e = np.power(y[:size] - u, 2) tmp = e[:int(size/avg_number)*avg_number] tmp.shape = -1, avg_number avg = np.average( tmp, axis=1 ) pl.plot(np.linspace(0, size, len(avg)), 10*np.log10(avg), linewidth=2.0, label=label) def diff_db(h0, h): return 10*np.log10(np.sum((h0-h)*(h0-h)) / np.sum(h0*h0)) # 用NLMS進行系統辨識的模擬, 未知系統的傳遞函數為h0, 使用的參照信號為x def sim_system_identify(nlms, x, h0, step_size, noise_scale): y = np.convolve(x, h0) d = y + np.random.standard_normal(len(y)) * noise_scale # 添加白色噪聲的外部干擾 h = np.zeros(len(h0), np.float64) # 自適應濾波器的長度和未知系統長度相同，初始值為0 u = nlms( x, d, h, step_size ) return y, u, h def system_identify_test1(): h0 = make_path(32, 256) # 隨機產生一個未知系統的傳遞函數 x = np.random.standard_normal(10000) # 參照信號為白噪聲 y, u, h = sim_system_identify(nlms_numpy.nlms, x, h0, 0.5, 0.1) print diff_db(h0, h) pl.figure( figsize=(8, 6) ) pl.subplot(211) pl.subplots_adjust(hspace=0.4) pl.plot(h0, c="r") pl.plot(h, c="b") pl.title(u"未知系統和收斂后的濾波器的系數比較") pl.subplot(212) plot_converge(y, u) pl.title(u"自適應濾波器收斂特性") pl.xlabel("Iterations (samples)") pl.ylabel("Converge Level (dB)") pl.show() def system_identify_test2(): h0 = make_path(32, 256) # 隨機產生一個未知系統的傳遞函數 x = np.random.standard_normal(20000) # 參照信號為白噪聲 pl.figure(figsize=(8,4)) for step_size in np.arange(0.1, 1.0, 0.2): y, u, h = sim_system_identify(nlms_numpy.nlms, x, h0, step_size, 0.1) plot_converge(y, u, label=u"μ=%s" % step_size) pl.title(u"更新系數和收斂特性的關系") pl.xlabel("Iterations (samples)") pl.ylabel("Converge Level (dB)") pl.legend() pl.show() def system_identify_test3(): h0 = make_path(32, 256) # 隨機產生一個未知系統的傳遞函數 x = np.random.standard_normal(20000) # 參照信號為白噪聲 pl.figure(figsize=(8,4)) for noise_scale in [0.05, 0.1, 0.2, 0.4, 0.8]: y, u, h = sim_system_identify(nlms_numpy.nlms, x, h0, 0.5, noise_scale) plot_converge(y, u, label=u"noise=%s" % noise_scale) pl.title(u"外部干擾和收斂特性的關系") pl.xlabel("Iterations (samples)") pl.ylabel("Converge Level (dB)") pl.legend() pl.show() def sim_signal_equation(nlms, x, h0, D, step_size, noise_scale): d = x[:-D] x = x[D:] y = np.convolve(x, h0)[:len(x)] h = np.zeros(2*len(h0)+2*D, np.float64) y += np.random.standard_normal(len(y)) * noise_scale u = nlms(y, d, h, step_size) return h def signal_equation_test1(): h0 = make_path(5, 64) D = 128 length = 20000 data = np.random.standard_normal(length+D) h = sim_signal_equation(nlms_numpy.nlms, data, h0, D, 0.5, 0.1) pl.figure(figsize=(8,4)) pl.plot(h0, label=u"未知系統") pl.plot(h, label=u"自適應濾波器") pl.plot(np.convolve(h0, h), label=u"二者卷積") pl.title(u"信號均衡演示") pl.legend() w0, H0 = scipy.signal.freqz(h0, worN = 1000) w, H = scipy.signal.freqz(h, worN = 1000) pl.figure(figsize=(8,4)) pl.plot(w0, 20*np.log10(np.abs(H0)), w, 20*np.log10(np.abs(H))) pl.title(u"未知系統和自適應濾波器的振幅特性") pl.xlabel(u"圓頻率") pl.ylabel(u"振幅(dB)") pl.show() signal_equation_test1() ```