dsp121.m

使用维纳滤波器对信号进行处理

p=input('请输入AR模型的阶数p=');
a=zeros(1,p);
for i = 1:p,                        %根据所输入的阶数输入参数
    teStr=sprintf('请输入AR模型的参数a%d=',i);
    a(i)=input(teStr);
end
tempA=zeros(1,p+1);                 %生成最高阶数为p的多项式，系数从高次
                                    %到低次依次为ap到a1，常数项为1
tempA(p+1)=1;
for i=1:p,
    tempA(p+1-i)=a(i);
end  
rootp=zeros(p,1);
rootp=roots(tempA)                   %求多项式的所有根
checkA=zeros(p,1);
for i=1:p,
checkA(i,1)=abs(rootp(i,1));       %求多项式的所有根的模
end
checkA
for i=1:p,             %检验多项式的所有根的模是否大于1，若全部大于1，继续执行下面的程序，否则退出
 if check(i,1)<=1; 
    display('s(n)不收敛，错误的输入,请重新输入');
    return;
    end
end
%p=input('请输入AR模型的阶数p=');
% a=zeros(1,p);
% for i = 1:p,
%     teStr=sprintf('请输入AR模型的参数a%d=',i);
%     a(i)=input(teStr);
% end
deltaW=input('请输入白噪声方差:');
L=input('请输入信号s(n)样本个数L:');
w=sqrt(deltaW)*randn(1,L);
s=zeros(1,L);
for i = 1:L,                                  %计算出s(n)
    te=0;
    for m = 1:p,
        if i-m>0,           
           te=te+a(m)*s(i-m);
        end
    end
    s(i)=w(i)-te;
end
Rss = zeros(p+1,p+1);                   %Rss为s(n)的自相关
for i = 1:p+1,                        %i为行数
    for j = 1:p+1,                    %j为列数
        teR=0;   
        te=abs(i-j);                %括号中的值是行数减去列数
        for k = 1:L-te,             %由式1-18得到Rss 
            teR=teR+s(k)*s(k+te);
        end
        teR=teR/(L-te);
        Rss(i,j)=teR;               %赋值给Rss
    end
end
A=zeros(p+1,1);
epsilon=zeros(p+1,1);
epsilon(1)=deltaW;
A=Rss^(-1)*epsilon;
for i = 1:p,     
    sprintf('a%d=%f\n',i,A(i+1))
end
deltaW_e=0;
for i = 1:p,
    deltaW_e=deltaW_e+a(i)*Rss(1,i+1);
end
deltaW_e=deltaW_e+Rss(1,1);
sprintf('deltaW_e=%f',deltaW_e)