南京邮电大学现代信号处理（MSP）课程大作业（2021）

之前在网上找的答案，要不就要下载券，要不就是PDF格式的，反正一点也不方便。所以我把可以成功运行的代码发到这里来，方便大家借鉴。

一共两个题目，一共是五个小题，老师的要求是五选二，所以这里第一题和第二题分别选做一小题。

第一题题目：

1. 根据所提供的数据(见群文件test.txt，该文件每相邻两个数组成一个复数据，一共可以组成64个复数据，它们表示了一个随机信号的64个取样值)，试用如下方法估计其功率谱，并画出不同阶数（10阶、25阶、40阶）情况下的功率谱曲线：

1. Levinson算法
2. Burg算法
3. MUSIC算法

解答：选做Levinson算法

原理啥的网上都可以找到，所以这里就放代码

function S=Levinson(X,p)
N=length(X) ;
for m=0:N-1
   R(m+1)=sum(conj(X([1 :N-m])).*X([m+1:N]))/N;
end
a=-R(2) /R(1) ;
sigma2=(1-abs(a)^2)*R(1) ;
gamma=-a;
for k=2:p
   sigma2(k)=R(1)+sum(a.*conj(R([2:k])));
   D=R(k+1)+sum(a.*R([k:-1:2]));
   gamma(k)=D/sigma2(k) ;
   sigma2(k)=(1-abs(gamma(k))^2)*sigma2(k);
   a=[a-gamma(k)*conj(fliplr(a)),-gamma(k)] ;
end
   sigma2=real(sigma2) ;
   f=linspace(-0.5,0.5,512) ;
for k=1:512
   S(k)=10*log10(sigma2(end)/(abs(1+sum(a.*exp(-j*2*pi*f(k)*[1:p])))^2));
end
plot(f,S) ;
title(['Levinson:',int2str(p),'阶']);
xlabel ('归一化频率'), ylabel ('功率谱/dB') ;

上面是函数体，放在m文件中

下面就是测试代码段，放在运行窗口就可以

load test.txt;
x=zeros(64,1);
for n=1:64
    x(n,1)=test(n,1)+j*test(n,2);
end
p=[10 25 40];
for k=1:3
   subplot(2,2,k);
   Levinson(x,p(k)) ;
end

当然test文件要放在同一个目录下，不然加载不了，test文件老师会给的，其实就是64行，每一行两个实数，要连接成复数，也就是64个复数。

结果大概是这样：

 第二题题目：
 

2. 图1为一个自适应均衡器的应用。其中，均衡器采用横向FIR结构（滤波器长度M=11）；系统输入是取值为±1的随机序列 ，其均值为零；参考信号 ；信道具有脉冲响应：

式中 用来控制信道的幅度失真(W = 2～4, 如取W = 2.9,3.1,3.3,3.5等)；信道受到均值为零、方差为 的高斯白噪声 的干扰（ 取值为0.01或者0.001时，相当于信噪比为20dB或30dB)。

   试比较基于下列两种算法的自适应均衡器在不同信道失真、不同的信噪比下的收敛情况(画出图1中误差 的均方值的学习曲线):

1. LMS算法
2. RLS算法

并分析其结果。

图1  自适应均衡器

解答：选做LMS算法

代码网上翻的

clear;
clf;
w=[2.9,3.1,3.3,3.5];
M=11;
sgm2=0.001; %噪声方差
T=7;
N=400; %训练次数
iteration=500;
L=iteration+T+M-1; %单个输入信号长度
u=0.035; %迭代步长
value=zeros(length(w),L-M+1-T) ;
for ww=1:length(w)
    h=ones(1,3);
    h(1)=0.5*(1+cos(2*pi/w(ww))) ;
    h(3)=h(1) ;
    e2=zeros(N,L-M+1-T) ;
    for n=1:N
        rand('seed',n*N) ;
        X=sign(2*rand(1,L)-1);%产生一个随机信号序列X=randsrc(1,L,[+1,-1], 1)
        D=zeros(size(X)) ;
        for mm=T+1:L 
            D(mm)=X(mm-T) ;
        end%对应的参考信号
        U=conv(X,h);
        randn('seed',n*N) ;
        V=randn(size(U)).*sqrt(sgm2); %产生高斯噪声
        R=U+V;%自适应滤波器输入信号
        W=zeros(M,1); %滤波器参数的初始值为0
        for m=1:iteration
           r=R(T+m:T+m+M-1)';
           y=r'*W;
           e=D(T+m+M-1)-y; %误差信号
           e2(n,m)=e.^2;
           W=W+2*u*e*r; %滤波器参数迭代
        end
    end
value(ww,:)=mean(e2); %均方误差值
end
semilogy(value(1,:),'k-') ;hold on;
semilogy(value(2,:),'b-.');
semilogy(value(3,:),'m:');
semilogy(value(4,:),'r--');
grid;
xlabel('重复次数');
ylabel('均方误差');
legend('W=2.9','W=3.1','W=3.3','W=3.5');

运行结果大概这样：