求助一个滤波器生成函数

峰雪无痕

　　在下面附件1中有待求的公式。
　　附件2是我写的程序，但是在有些数据时会出错
　　比如：
　　K=16;
　　M=190;
　　B=10;
　　t=0:1/(M*K):(K-1); %序列横座标
　　N=length(t); %序列总长
　　fs=1/(M*K);
　　n1=randn(1,N);
　　t1=t;
　　b=[0.2,0.4,0.6,0.8,1.0,1.0,1.0,0.7,0.5,0.3,0.1];
　　[n2,tn2]=xgzs(n1,t1,b,fs);
　　n2=B*n2;

峰雪无痕

因为可能存在看不到附件的情况，我把程序在下面写出来，请大家帮忙看下！

function [n2,tn2]=xgzs(n1,tn1,b,fs)
% 知道噪声n1(t)的情况下求其相关噪声n2(t)
zc=length(b);                           %序列b(n)的长度
n2=0;
t2=0;
for i=0:1:(zc-1)
        n0=n2;
    t0=t2;
    [n3,t3]=Delay_sig0(n1,tn1,i,fs);
        n2=b(i+1)*n3;
    t2=t3;   
        [n2,t2]=Sig_Add(n0,t0,n2,t2,1,1,fs);
    i
end
n2=n2;
tn2=t2;

峰雪无痕

function [y,n]=Delay_sig(x,m,n0,fs)                %  n0 为延迟的样点数
n=m-n0*fs;
y=x;

function [y,n]=Sig_Add(x1,n1,x2,n2,a,b,fs)
n=min(min(n1),min(n2)):fs:max(max(n1),max(n2));
y1=zeros(1,length(n));
y2=y1;
y1(find((n>=min(n1))&(n<=max(n1))==1))=x1;
y2(find((n>=min(n2))&(n<=max(n2))==1))=x2;
y=a*y1+b*y2;