ADSP现代数字信号处理仿真实验报告

1、目 录仿真一：LMS算法和RLS算法11 自适应滤波的基本原理11.1 自适应最小均方（LMS）算法11.2 递归最小二乘方（RLS）算法22 仿真实验43 结果分析6仿真二：P阶Levinson-Durbin算法81 要求：82 算法描述82.1 产生信号82.2 L-D算法92.3 对比信号谱功率和LD算法谱估计103 结果分析113.1 AR模型113.2 MA模型123.3 总结13仿真一：LMS算法和RLS算法1 自适应滤波的基本原理自适应滤波器由参数可调的数字滤波器/自适应处理器和自适应算法两部分组成，如图1所示。输入信号x(n)通过参数可调数字滤波器后产生的输出信号为y(n)，将

2、其与参考信号d(n)进行比较，得到误差信号e(n)。误差信号e(n)经过一定的自适应算法后反馈到参数可调数字滤波器，对滤波器进行参数调整（有时还需要利用x(n)），以使得e(n)最终的均方值最小。这是一种自动控制理论，因此，滤波器在设计时不需要事先知道输入信号和噪声的统计特性，而能够根据输入信号的统计特性变化自动跟踪这种变化，自动调整参数，使滤波器性能达到最佳。图 1 自适应滤波器框图图1所示自适应滤波器，输入信号为：x(n)和d(n)，两个输出为：y(n)和e(n)。当误差信号e(n)的均方误差达到最小的时候，可以证明信号y(n)是信号d(n)的最佳估计。1.1 自适应最小均方（LMS）算法

3、最陡下降法每次迭代都需要知道性能曲面上某点的梯度值，而梯度值只能根据观测数据进行估计。LMS算法是一种有用简单的估计梯度的方法，其最核心的思想是采用平方误差最小代替均方误差最小准则。信号基本关系：式中，W(n) 为 n 时刻自适应滤波器的权矢量，下一时刻权矢量 W(n+1) 等于当前权矢量 Wn 加上一个修正量，该修正量是误差信号en的加权值，加权系数 2x(n) 正比于当前的输入信号 x(n)。 N 为自适应滤波器的阶数；x(n) 为 n时刻自适应滤波器的参考输入矢量，由最近N 个信号采样值构成，；d(n)是期望的输出值；e(n) 为自适应滤波器的输出误差调节信号； 是控制搜索步长的参数称为

4、自适应增益常数，又叫收敛因子或步长因子。算法实现如下：步骤一：信号生成（点数取为500），其matlab代码如下：clcclose all;clear all; n=500; %设置信号点数a1=-1.6;a2=0.8; %模型参数赋值k=randn(1,n)' %生成白噪声信号x=zeros(1,n)'x(1)=k(1);x(2)=k(2)-a1*x(1);for i=3:n %生成信号x(n) x(i)=k(i)-a1*x(i-1)-a2*x(i-2);end;plot(x,'b-');%画出信号和噪声波形hold on;plot(k,'r-'

5、;);% LMS算法 %L=2;%自回归模型长度u=0.002;%LMS自适应参数取值0.002w=zeros(L,n);%LMS滤波器的系数for i=(L+1):n X=x(i-1:-1:(i-L); y(i)=X'*w(:,i); e(i)=x(i)-y(i); w(:,(i+1)=w(:,i)+2*u*e(i)*X;end;a1=-w(1,:);m=1:n;figure();plot(m,a1(m),'b-');hold on;1.2 递归最小二乘方（RLS）算法递归最小二乘方（RLS）算法是FIR维纳滤波器的一种时间递归算法，严格以最小二乘方准则为依据，收敛速

6、度快。RLS算法的关键是用二乘方的时间平均的最小化准则取代最小均方准则，并按时间进行迭代计算。其基本原理如下所述，其中称为遗忘因子，它是小于等于1的正数；为参考信号或期望信号；为第n次迭代的权值；均方误差。首先，对初始时刻到当前时刻所有误差的平方进行平均并使其最小化，然后在这一准则指导下确定FIR滤波器的权系数矢量即越旧的数据对的影响越小。对滤波器权系数矢量求偏导数，并令结果等于零整理得到标准方程定义标准方程可以化简成形式：经求解可以得到迭代形式定义：，则可知T的迭代方程为系数的迭代方程为其中增益和误差的定义分别为总结RLS算法步骤如下：1) 在时刻n，已经知道和也已经存储在滤波器的延时部件中

7、。2) 利用公式（1.9）、（1.10）、（1.11）和（1.12）计算，并得到滤波器的输出响应和误差，即：3) 进入第次迭代。算法实现如下：% RLS滤波 %L=2;lamda=1;%赋值，初值为1w=zeros(L,n);%权系数，初值为0T=eye(L,L)*10;%T初始值为10for i=(L+1):n X=x(i-1:-1:(i-L); K=(T*X)/(lamda+X'*T*X); e1=x(i)-w(:,i-1)'*X; w(:,i)=w(:,i-1)+K*e1; y(i)=w(:,i)'*X; e(i)=x(i)-y(i); T=(T-K*X'

8、*T)/lamda;end;a1=-w(1,:);plot(m,a1(m),'r');hold on;convergence=-1.6*ones(1,n);%系数收敛到-1.6plot(m, convergence(m),'g:');%绘制收敛参考线axis(1,n,-2,0);grid on;title('自适应权系数a1(n)的过渡过程（LMS & RLS 算法比较）')xlabel('n');ylabel('a1(n)');legend('LMS算法','RLS算法')

9、;2 仿真实验仿真过程按照如下过程进行1) 信号生成。首先产生零均值、单位方差的高斯白噪声序列，然后将此通过一个简单的二阶自回归滤波器生成信号。二阶自回归滤波器的模型，参数设置为，。2) 将上一步骤生成的信号通过LMS和RLS自适应滤波器进行处理。3) 绘制各种图形曲线。仿真结果如下：1) 信号和高斯白噪声的波形如图 2所示：图 2 数据和高斯白噪声波形2) 将信号输入LMS和RLS滤波器进行处理，迭代之后生成的图形，如图 3所示。图 3 自适应权系数的过渡过程（LMS算法和RLS算法收敛曲线对比）3) 最后画出的取值对RLS算法的影响。代码如下：figure();convergence=-1

10、.6*ones(1,n); plot(m, convergence(m),'g:'); axis(1,n,-2,0);title('遗忘因子对RLS算法的影响');xlabel('n');ylabel('a1(n)');L=2;lamda=0.98; %赋值为0.98w=zeros(L,n); T=eye(L,L)*10; for i=(L+1):n X=x(i-1:-1:(i-L); K=(T*X)/(lamda+X'*T*X); e1=x(i)-w(:,i-1)'*X; w(:,i)=w(:,i-1)+K*e1

11、; y(i)=w(:,i)'*X; e(i)=x(i)-y(i); T=(T-K*X'*T)/lamda;end;a1=-w(1,:);plot(m,a1(m),'g');遗忘因子时，RLS收敛速度变慢，且收敛性变差，如图4所示。图 4 遗忘因子对RLS收敛速度的影响（）3 结果分析通过实验仿真发现，RLS算法收敛速度明显优于LMS算法，如图3所示。RLS算法系数在的时候可以很快的收敛到最佳值，而LMS算法则需要花费更多的时间。取值越小，RLS算法的收敛速度越慢，这是因为现在有效记忆长度只有，在计算最佳权系数时没有充分利用能够获得的全部取样数据。仿真二：P阶Le

12、vinson-Durbin算法1 要求：A. 使用AR(2)模型（b0 = 1, a1 = 0, a2 = 0.81）和MA(2)模型（bn = 1, 1, 1）处理1000个样本点B. 对AR(2)和MA(2)使用L-D求解模型参数，p分别为2和10C. 画出这两种情况下产生的AR谱，并列出反射系数表D. 描述L_D算法，计算一个P阶L_D算法的计算量2 算法描述2.1 产生信号输入激励u(n)是均值为零、方差为2的白噪声序列，其输出与输入之间满足差分方程xn= -k=1pakxn-k+k=0qbku(n-k) 其中bk是前馈支路的系数，称为MA系数；ak是反馈支路的系数，称为AR系数。由仿

13、真要求，AR(2)模型参数为，b0 = 1，a1 = 0，a2 = 0.81，则模型输入输出方程为xn= -0.81xn-2+u(n)MA(2)模型参数为bn= 1,1,1，则模型输入输出方程为xn=un+un-1+u(n-2)产生两个信号的代码如下所示：% signal productionN = 1000; % sampling numberu = randn(N,1); % white noise%AR(2)a = 0, 0.81;AR = zeros(N,1);AR(1) = u(1);AR(2) = u(2) - a(1) * AR(1);for i = 3: N AR(i) = -

14、a(1) * AR(i-1) -a(2) * AR(i-2) + u(i);end%MR(2)b = 1, 1, 1;MA = zeros(N,1);MA(1) = b(1) * u(1);MA(2) = b(1) * u(2) + b(2) * u(1);for i = 3: N MA(i) = b(1) * u(i) + b(2) * u(i-1) + b(3) * u(i-2);end2.2 L-D算法L_D算法是一种按阶次进行递推的算法，即首先以AR(0)和AR(1)模型参数作为初始条件，计算AR(2)模型参数，再以此计算AR(3)的参数，一直到计算出AR(p)的参数为止。推导结果如下

15、：k2=i=0kak,iRi,ak,0=1Dk=i=0kak,iR(k+1-i),ak,0=1k+1=Dkk2k+12=k2-k+1Dkak+1,i=ak,i-k+1ak,k+1-i第一、二个公式需要的计算次数为2k+1次乘法和2k次加法；第三个公式需要一次除法；第四个公式需要一次加法一次乘法；第五个公式需要k次乘法和k次加法；因而由k阶计算k+1阶系数总共需要3k+2次乘法，3k+1次加法和一次除法；因而可知，计算p阶的系数总共所需的乘法数目为： k=0p-1(3k+2)=3p2+P2次乘法， k=0p-1(3k+1)=3p2-P2次加法和p次除法。通过以上结果公式，设计如下LD函数，输入为

16、自相关系数R，输出AR/MA系数，激励源标准差和反射系数。function a_k, sigma_k, gama_k= LD ( R )% input: R-coefficient of autocorrelation% output: a_k-coefficient of AR model% sigma_k-standard deviation% gama_k-reflection coefficientp = length(R) - 1;a_k = zeros(p+1, 1);sigma_k = zeros(p+1, 1);gama_k = zeros(p+1, 1); sigma_k(1

17、) = R(1);D(1) = R(2);gama_k(1) = D / sigma_k(1);a_k(1) = 1;a_k(2) = -gama_k(1); for i = 1:p-1 sigma_k(i+1) = R(1)+a_k(2:(i+1)'*R(2:(i+1); D(i+1) = a_k(1:(i+1)'*R(i+2):-1:2); gama_k(i+2) = D(i+1)/sigma_k(i+1); sigma_k(i+2) = (1-gama_k(i+2)2)*sigma_k(i+1); a_k(i+2) = -gama_k(i+2); a_k(2:(i+1)

18、= a_k(2:(i+1)-gama_k(i+1)*a_k(i+1):-1:2);end end2.3 对比信号谱功率和LD算法谱估计通过产生信号，计算信号功率谱，再通过LD算法计算AR模型参数，估算功率谱，并将二者进行比较，从而得出谱估计精度。代码如下：% L-D algorithm% AR(2) parameter: P = 2p = 2;k = 0:2;w = 0:0.01:2*pi;A = 1 0 0.81'H_AR = Sigma2 ./(abs(A'*exp(-j*k'*w).2);% power spectrum % R(m) = Ex(n)x(m+n)k

19、 = 0:p;for i = 0: p R_AR(i+1) = AR(1:N-i)' * AR(i+1:N) / N;enda_AR, sigma_AR, gama_AR = LD( R_AR' );H_LD1 = sigma_AR(p+1)./ (abs(A'*exp(-j*k'*w).2); figure(1)plot(w, H_AR, 'r-', w, H_LD1, 'b-');xlabel('w');ylabel('|H_AR|');title('AR模型功率谱'); % MA(2) parameter: P = 10q = 10;k=0:2;B=1 1 1'H_MA=abs(B'*exp(-j*k'*w).2; % R(m) = Ex(n)x(m+n)for i=0:q R_MA(i+1) = MA(1:N-i)'*MA(i+1:N)/N;endk