LMS算法仿真（Matlab）

william

1. 程序1：基本LMS算法
   
2. %  该程序实现时域LMS算法，并用统计的方法仿真得出不同步长下的收敛曲线
   
3. clear                             %   清空变量空间
   
4. g=100;                          %   统计仿真次数为g
   
5. N=1024;                        %   输入信号抽样点数N
   
6. k=128;                          %   时域抽头LMS算法滤波器阶数
   
7. pp=zeros(g,N-k);                 %   将每次独立循环的误差结果存于矩阵pp中，以便后
   
8.                                %   面对其平均
   
9. u=0.001;                        %   收敛因子
   
10. for q=1:g               
    
11. t=1:N;
    
12.     a=1;
    
13.     s=a*sin(0.5*pi*t);            %   输入单频信号s
    
14.     figure(1);
    
15.     plot(t,real(s));                %   信号s时域波形
    
16.     title('信号s时域波形');
    
17.     xlabel('n');
    
18.     ylabel('s');
    
19.     axis([0,N,-a-1,a+1]);
    
20.     xn=awgn(s,3);               %   加入均值为零的高斯白噪声，信噪比为3dB   
    
21.     % 设置初值
    
22.     y=zeros(1,N);                %    输出信号y
    
23.     y(1:k)=xn(1:k);              %   将输入信号xn的前k个值作为输出y的前k个值
    
24.     w=zeros(1,k);                %    设置抽头加权初值
    
25.     e=zeros(1,N);                %    误差信号
    
26.     % 用LMS算法迭代滤波
    
27.     for i=(k+1):N
    
28.                XN=xn((i-k+1):(i));
    
29.                y(i)=w*XN';
    
30.                e(i)=s(i)-y(i);
    
31.                w=w+u*e(i)*XN;   
    
32.      end
    
33.      pp(q,:)=(e(k+1:N)).^2;  
    
34.   end
    
35.   for b=1:N-k
    
36.          bi(b)=sum(pp(:,b))/g;     %   求误差的统计平均
    
37.   end
    
38.   figure(2);                    %   算法收敛曲线
    
39.   t=1:N-k;
    
40.   plot(t,bi);
    
41.   hold off                     %   将每次循环的图形显示结果保存下来
    
42. 
    
43. 程序2  归一化LMS算法
    
44. 
    
45. MATLAB程序实现如下：
    
46. 1. NLMS算法1次实验
    
47. %  N=训练序列长度
    
48. %  u=收敛因子
    
49. clear;
    
50. N=500;
    
51. db=20;
    
52. sh1=sqrt(10^(-db/10));
    
53. u=1;
    
54. error_s=zeros(1,N);
    
55. for  loop=1:1
    
56.         w=0.05*ones(1,11)';
    
57.         V=sh1*randn(1,N );
    
58.         K=randn(1,N)-0.5;
    
59.         x=sign(K);
    
60. for n=3:N;
    
61.             M(n)=0.3*x(n)+0.9*x(n-1)+0.3*x(n-2);
    
62. end
    
63. z=M+V;
    
64. for n=8:N;
    
65.             d(n)=x(n-7);
    
66. end
    
67.         a(1)=z(1)^2;
    
68. for n=2:11;
    
69.             a(n)=z(n).^2+a(n-1);
    
70. end
    
71. for n=12:N;
    
72.             a(n)=z(n).^2-z(n-11)^2+a(n-1);
    
73. end
    
74. for n=11:N;
    
75.             z1=[z(n) z(n-1) z(n-2) z(n-3) z(n-4) z(n-5) z(n-6) z(n-7) z(n-8) z(n-9) z(n-10)]';
    
76.             y(n)=w'*z1;
    
77.             e(n)=d(n)-y(n);                                                                           
    
78.             w=w+u./(eps+a(n)).*z1.*conj(e(n));
    
79. end                                                                                       
    
80.         error_s=error_s +e.^2;
    
81. end
    
82. w
    
83. error_s=error_s./1;
    
84. n=1:N;
    
85. plot(n,error_s);
    
86. xlabel('n (当u=1;DB=20时)');                                                                              
    
87. ylabel('e(n)^2');
    
88. title('NLMS算法1次实验误差平方的均值曲线');
    
89. 2.NLMS算法20次实验
    
90. clear;
    
91. N=500;
    
92. db=20;
    
93. sh1=sqrt(10^(-db/10));
    
94. u=1;
    
95. error_s=zeros(1,N);
    
96. for   loop=1:20
    
97.         w=0.05*ones(1,11)';
    
98.         V=sh1*randn(1,N );
    
99.         K=randn(1,N)-0.5;
    
100.         x=sign(K);
     
101. for  n=3:N;
     
102.            M(n)=0.3*x(n)+0.9*x(n-1)+0.3*x(n-2);
     
103. end
     
104.         z=M+V;
     
105. for n=8:N;
     
106.              d(n)=x(n-7);
     
107. end
     
108.         a(1)=z(1)^2;
     
109. for n=2:11;
     
110.              a(n)=z(n).^2+a(n-1);
     
111. end
     
112. for n=12:N;
     
113.             a(n)=z(n).^2-z(n-11)^2+a(n-1);
     
114. end
     
115. for n=11:N;
     
116.             z1=[z(n) z(n-1) z(n-2) z(n-3) z(n-4) z(n-5) z(n-6) z(n-7) z(n-8) z(n-9) z(n-10)]';
     
117.             y(n)=w'*z1;
     
118.             e(n)=d(n)-y(n);                                                                           
     
119.             w=w+u./(eps+a(n)).*z1.*conj(e(n));
     
120. end                                                                                       
     
121.         error_s=error_s +e.^2;
     
122. end
     
123. w
     
124. error_s=error_s./20;
     
125. n=1:N;
     
126. plot(n,error_s);
     
127. xlabel('n (当u=1;DB=20时)');                                                                              
     
128. ylabel('e(n)^2');
     
129. title('NLMS算法20次实验误差平方的均值曲线');
     
130. 
     
131. 
     
132. 
     
133. 基于LMS算法的系统辨识
     
134. 
     
135. clear
     
136. clc
     
137. ee=0;
     
138. fs=800;
     
139. det=1/fs;
     
140. f1=100;
     
141. f2=200;
     
142. t=0:det:2-det;
     
143. x=randn(size(t))+cos(2*pi*f1*t)+cos(2*pi*f2*t);
     
144. %未知系统
     
145. [b,a]=butter(5,150*2/fs);
     
146. d=filter(b,a,x);
     
147. %自适应FIR滤波器
     
148. N=5;
     
149. delta=0.06;
     
150. M=length(x);
     
151. y=zeros(1,M);
     
152. h=zeros(1,N);
     
153. for  n=N:M
     
154.          x1=x(n:-1:n-N+1);
     
155.          y(n)=h*x1';
     
156.          e(n)=d(n)-y(n);
     
157.          h=h+delta.*e(n).*x1;
     
158. end
     
159. X=abs(fft(x,2048));
     
160. Nx=length(x);
     
161. kx=0:800/Nx:(Nx/2-1)*(800/Nx);
     
162. D=abs(fft(d,2048));
     
163. Nd=length(D);
     
164. kd=0:800/Nd:(Nd/2-1)*(800/Nd);
     
165. Y=abs(fft(y,2048));
     
166. Ny=length(Y);
     
167. ky=0:800/Ny:(Ny/2-1)*(800/Ny);
     
168. figure(1);
     
169. subplot(3,1,1)
     
170. plot(kx,X(1:Nx/2));xlabel('Hz')
     
171. title('原始信号频谱')
     
172. subplot(3,1,2)
     
173. plot(kd,D(1:Nd/2))
     
174. title('经未知系统后信号频谱');xlabel('Hz')
     
175. subplot(3,1,3)
     
176. plot(ky,Y(1:Ny/2))
     
177. title('经自适应FIR滤波器后信号频谱');xlabel('Hz')

复制代码

转自：http://blog.sina.com.cn/s/blog_5def5a660100c25e.html

sizhiyuan2006

谢谢分享！！！