wholework.m

内插滤波器

pntaps=[0 0 1 0 0 0 0 0 0 1];
pninitial=[0 0 0 0 0 0 0 1 0 1];
pndata=zeros(1,1023);
pnregister=pninitial;
n=0;
k=0;
while k==0
    n=n+1;
    pndata(1,n)=pnregister(1,1);
    feedback=rem((pnregister*pntaps'),2);
    pnregister=[feedback,pnregister(1,1:9)];
    if pnregister==pninitial; k=1; end
end
pndata=2*pndata-1;

hbfnum=[2.2573684192062548e-003,-2.9949007724012425e-018,-6.5559151107721057e-003,6.0279437638757836e-018,1.4128851170508742e-002,-9.6339260869033891e-018,-2.6841334662664040e-002,1.3327489538311222e-017,4.9085258153930590e-002,-1.6235012868626662e-017,-9.6807091322389766e-002,1.8717245638673309e-017,3.1511491786405837e-001,5e-001,3.1511491786405837e-001,1.8717245638673309e-017,-9.6807091322389766e-002,-1.6235012868626662e-017,4.9085258153930590e-002,1.3327489538311222e-017,-2.6841334662664040e-002,-9.6339260869033891e-018,1.4128851170508742e-002,6.0279437638757836e-018,-6.5559151107721057e-003,-2.9949007724012425e-018,2.2573684192062548e-003];
hbf1=zeros(1,1023*2);
m=0
while m<1023
    m=m+1;
    hbf1(1,(m-1)*2+1)=pndata(1,m);
end
hbf_1=conv(hbf1,hbfnum);

hbf2=zeros(1,2072*2);
m=0
while m<2072
    m=m+1;
    hbf2(1,(m-1)*2+1)=hbf_1(1,m);
end
hbf_2=conv(hbf2,hbfnum);

hbf3=zeros(1,4170*2);
m=0
while m<4170
    m=m+1;
    hbf3(1,(m-1)*2+1)=hbf_2(1,m);
end
hbf_3=conv(hbf3,hbfnum);


k=0;
m=0;
while m<8366
    m=m+1;
    if m<=2; cicdata(1,m)=hbf_3(1,m); 
    else
    cicdata(1,m)=hbf_3(1,m)-hbf_3(1,m-2);
end
end
cicin=zeros(1,8366*7);
k=0;
m=0;
while m<8366
    m=m+1;
    cicin(1,(m-1)*7+1)=cicdata(1,m);
end

cicout=zeros(1,8366*7);
m=0;
while m<8366*7
    m=m+1;
    if m<2; cicout(1,m)=cicin(1,m);
    else 
    cicout(1,m)=cicin(1,m)+cicout(1,m-1);
end
end



i=[1:1:8366*7];
stem(i,cicout);axis([1 8366*7 -10 10]);

n=[1:8366*7];
k=[-200:200];
w=(pi/100)*k;
x=cicout;

X=x*(exp(-j*pi/100)).^(n'*k);
magX=abs(X);angX=angle(X);
subplot(2,1,1);plot(w/pi,magX,'r');grid;
xlabel('\omega(角频率)/\pi');ylabel('|X|');title('幅度');
subplot(2,1,2);plot(w/pi,angX/pi,'r');grid;
xlabel('\omega(角频率)/\pi');ylabel('\times\pi rads');title('相位');