关于重采样的问题，麻烦大神帮忙看看

自己在matlab里面编的一个重采样程序，但是效果很差，而且输出前一段数据有错，麻烦路过的大神帮我看看？

%input data
; L9 s( d- H/ V6 s' C4 D# cfa = 8000; %%signal frequency
# a, u. K  S* F0 nfs = 44100; %%44.1kHz sampling frequency1 K  k4 N2 L3 c6 u
n = 1:64;
5 x! t* T8 q* `9 Ux = sin(2*pi*fa*n/fs);
0 G/ Q0 m0 t( h: n2 H' Mlengthx = length(x);
: d; \  O; h) a6 ~4 a1 n9 wt = n*1000/128/fs;* E& Z. O& B+ k0 T4 ]. X- T
%plot(t,x);xlabel('time in ms');

% [y fs nbits] = wavread('acoustic.wav');
' H9 i/ ~5 [8 |6 K% x = y';: V- z- q/ b/ b  [3 b6 {
% lx = length(x);
; [4 r) b3 v4 F' [1 z8 G%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
" k8 {! K) K$ Z$ j0 N9 h' k. v; d%Resample parameters
  |. A. p* j, ]5 FL1 = 8;! }# z9 |9 U3 W
M1 = 7;
1 E! ~: X0 m4 {( F0 z" O! iL2 = 4;
0 \" e, L9 k+ u+ Z, P  c% }. T" cM2 = 3;
/ h) z; V% o% G& y' g. {9 yL3 = 10;
/ V9 t$ M4 u: v+ |M3 = 7;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
. m: W* F/ X* H- w1 }  d%%filters
. O7 L4 n8 \& ^) Y( x+ eload filter.mat
- Z* _4 L0 |3 B( r1 jfilter1 = b1;
' y7 T8 J* `; B, J) lfilter2 = b2;) w7 I& Y+ p9 S. v* e
filter3 = b3;
, d# f1 O  l# A! }+ J%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
& C+ ~! x4 V  J& K%Initialize the subfilters

%First Filter
- t) o6 e$ _: `) W# l7 V& cupFilter1 = decompose(filter1,L1);
$ H4 |, j  a. F$ qdownFilter1 = decompose(filter1,M1);

%Second Filter" ]4 {6 d8 Q6 x9 f
upFilter2 = decompose(filter2,L2);; Z3 w5 Z- v, T( o* V/ j
downFilter2 = decompose(filter2,M2);

%Second Filter" c* X9 g' ?8 q* K
upFilter3 = decompose(filter3,L3);
0 y. Q6 }: U! E  N* v/ `6 B- ~! ydownFilter3 = decompose(filter3,M3);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
9 o' z" x  o3 H9 s: W/ j3 D%Polyphase Filtering

%First Filter
7 x: \- |( o- J" v; {" B%
+ u+ G$ q4 b4 B* D; Q, N/ A%%upsampling
: R9 x8 ?+ O2 `; \/ L& ointermediatex1 = zeros(L1,length(x));! G: M0 l5 P4 X+ D) Z
x11 = zeros(L1,length(x)*L1);8 _! a. c* ?3 s' E3 p
for i = 11
9 \4 y3 U! C5 c- z# {intermediatex1(i, = filter(upFilter1(i,,1,x); %%do the filtering before upsampling8 ^  \/ y4 T# i9 R
x11(i, = upsample(intermediatex1(i,:),L1); %%upsample by L, inserting L-1 zeros.# f# N0 {- w& g7 W8 B* V( q; p9 n
end

X1 = zeros(1,length(x)*L1);* l1 u! \9 ]( ~* {6 T
for i = 1:length(x)
! c7 t( O* L5 ?/ K" hX1(L1*(i-1)+11*(i-1)+L1) = x11(:,(i-1)*L1+1)';1 A7 z6 u+ D7 g' o, n7 t" W
end
3 W. H: a0 y$ l%%downsampling
# ^7 w4 [3 e0 L& sX1_down = [X1 zeros(1,M1-mod(length(X1),M1))];: Y3 d, I9 d" j8 s3 f2 ]
output1 = zeros(1,length(X1_down)/M1);
% }2 v* b( e" }x1down = zeros(M1,length(X1_down)/M1);
. L/ w8 e" V7 i  R) K, Xx11down = x1down;7 [, M/ D& j  ^4 {# v
for i = 1:M1
2 V& e. u# C- d& I* ^) O( s5 t& xfor k = 1:length(X1_down)/M1
. E, ?0 F! c; J6 I7 u: p# H3 Sx1down(i,k) = X1_down(M1*(k-1)+i);+ P+ R6 r6 N# V6 ~! s  k" N- t
end
; B% R5 ?- e1 S; j9 i2 D4 N; g% i7 Bx11down(i,:) = filter(downFilter1(1,:),1,x1down(i,:));, }5 H6 k* a/ _1 {4 g
output1 = output1 + x11down(i,:);* D: p  ^- {0 }  }0 Q
end

%Second Filter4 n' s# x+ k1 R9 N( u' t
%
; F9 _  u' c& \5 V4 m0 ~! j0 n3 ~/ m- j%%upsampling& [# a9 I8 W4 Y. u4 C7 s$ v
intermediatex2 = zeros(L2,length(output1));
" ]3 J( I; \; s+ x4 jx22 = zeros(L2,length(output1)*L2);8 o4 U" k/ R) h; f1 S
for i = 121 [6 `# h% E; R- P# S
intermediatex2(i,:) = filter(upFilter2(i,:),1,output1); %%do the filtering before upsampling1 B5 H7 h1 j6 [
x22(i,:) = upsample(intermediatex2(i,:),L2); %%upsample by L, inserting L-1 zeros.
1 G8 T; ?1 F" h" w" Gend

X2 = zeros(1,length(output1)*L2);
2 o0 Y! @$ h! q: `/ a4 q  i+ kfor i = 1:length(output1)) a! d4 H+ ]2 }$ e+ m
X2(L2*(i-1)+1:L2*(i-1)+L2) = x22(:,(i-1)*L2+1)';% S+ G0 ~0 v' `
end9 N: f. G+ s4 [/ ^
%%downsampling
5 r- v' @! }$ P+ lX2_down = [X2 zeros(1,M2-mod(length(X2),M2))];
7 ?0 f3 Y1 ^$ d) o9 }& @output2 = zeros(1,length(X2_down)/M2);
6 D+ I/ T' o7 ?9 Cx2down = zeros(M2,length(X2_down)/M2);
: G0 A& y, b! l: }- z6 d- S5 R3 jx22down = x2down;
. Q, W3 F( m# T7 E$ ?6 dfor i = 1:M2! c1 z8 M! j  \7 V! K
for k = 1:length(X2_down)/M22 G' z8 i( |) }, S6 D
x2down(i,k) = X2_down(M2*(k-1)+i);! z; M7 ~: k" R4 x7 ^9 I
end
1 t0 U+ ^0 y+ e9 Rx22down(i,:) = filter(downFilter2(1,:),1,x2down(i,:));
0 W/ T; U) i* b7 `0 Ioutput2 = output2 + x22down(i,:);
- T% e5 v) v, N- Lend

%Third Filter% O2 a( P0 f6 t7 r
%! Q! ~, o' `6 T+ l* h' e
%%upsampling( Y% [* n) V( S
intermediatex3 = zeros(L3,length(output2));8 b4 s8 p% ?" A
x33 = zeros(L3,length(output2)*L3);
) L4 p0 W% l2 N7 nfor i = 1:L3
, F, o9 B% @) ?/ h  G! Nintermediatex3(i,:) = filter(upFilter3(i,:),1,output2); %%do the filtering before upsampling
5 L9 T% v: V5 M4 ]1 n; j5 |x33(i,:) = upsample(intermediatex3(i,:),L3); %%upsample by L, inserting L-1 zeros.
2 c! `$ f9 m8 q0 Y* b9 aend

X3 = zeros(1,length(output2)*L3);
6 E7 m) f$ |  N2 |  gfor i = 1:length(output2)
) G9 o% P7 @8 I: tX3(L3*(i-1)+1:L3*(i-1)+L3) = x33(:,(i-1)*L3+1)';
4 W: Y  |- ?- i# @  f9 [end
( j/ R5 k0 b" n%%downsampling; Q8 @8 b# H. u( }5 D+ S
X3_down = [X3 zeros(1,M3-mod(length(X3),M3))];
8 U7 a2 _. s+ D  p/ ?6 {& xoutput3 = zeros(1,length(X3_down)/M3);3 S" }4 M/ }  V) f. Z; ~) \! D
x3down = zeros(M3,length(X3_down)/M3);
" H2 f$ u# a! a; K( v4 y$ Ax33down = x3down;
, q# H, P, D8 J6 Xfor i = 1:M38 M8 o/ c. x, L3 t
for k = 1:length(X3_down)/M39 D( I: n4 e9 J6 S$ c3 E
x3down(i,k) = X3_down(M3*(k-1)+i);
1 J. m- L$ p2 l: q8 C& f, Lend1 x. |" ~+ |5 M, i! j/ E
x33down(i,:) = filter(downFilter3(1,:),1,x3down(i,:));
! r' K) T- U! O2 q% I* J: ooutput3 = output3 + x33down(i,:);0 E; V8 q8 T9 \+ Q) L, F; A
end

output3 = output3/max(output3);0 ~7 r5 G1 U/ h8 D) P
lout = length(output3);

nx = 1:length(output3);
9 y5 V' [+ v8 ]xorigin = sin(2*pi*fa*nx/96000);4 W' c. _8 w4 b; ^# U  h* S$ F
error1 = xorigin - output3;
/ \8 \. O3 b, J0 Q8 fyin = [0 0 resample(x,320,147)];4 z; g) b( L4 _; ^( U  ~& p
error2 = xorigin - yin;

figure(1)
. Q3 F$ H+ @- J/ c  x* vsubplot(3,1,1);plot(xorigin);title('True Signal');2 l) j( o. V4 f) p( H. R
subplot(3,1,2);plot(output3);title('Resampled Signal');) V2 {/ J- O0 g' O1 G
subplot(3,1,3);plot(error1);title('Error');

figure(2)
1 _7 L( C0 R! vsubplot(3,1,1);plot(xorigin);title('True Signal');8 D9 r7 u! ]$ u  b: e
subplot(3,1,2);plot(yin);title('Built-in Resampled Signal');
' f# s. O1 o9 }3 r: w0 @subplot(3,1,3);plot(error2);title('Error');

function x_decompose = decompose(x,factor)
$ P0 x! \6 D! w8 V  l* {lx = length(x);

if mod(lx,factor) ~= 0;
1 o* W9 J# f1 b+ a+ {( ^( ^6 X7 yX = [x zeros(1,factor-mod(lx,factor))];
  {5 Z* l* j4 }/ t5 welse
) z( p# }8 X* \  T& }/ yX = x;
# [1 l" `& I8 e; N( O8 Send

x_decompose = zeros(factor,length(X)/factor);; U4 Z4 t5 V* ~% b9 Q6 o
for i = 1:factor
( R/ b( w, ^* R# Vfor k = 1:length(X)/factor" R" N$ `+ b3 s9 @
x_decompose(i,k) = X(factor*(k-1)+i);' Q; S$ I2 n* k5 q
end& D. z" h2 K& I  a
end

end

%Filter Generation
' m: o' Y5 n- |) ]0 Sfs = 44100;
9 u2 t3 u# I' i5 ]L1 = 8;/ m- a  W+ M: P* k, s8 f5 t. V
M1 = 7;
& j( ^; y5 \& Y$ L* L+ E9 sL2 = 4;
. F3 O& ], t3 ^2 EM2 = 3;' n6 V- t7 ?2 z& X8 B- t
L3 = 10;2 Y. d* `: ]) i$ w3 ]  V
M3 = 7;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%6 n- {/ ~3 {# S9 v0 `( ?
%+ ^/ T- @& m& V! O- \5 H
%%First Filter' T4 F7 o+ X- U+ y# M9 w
f1 = [20000 24000];
5 b" z4 f( N8 D4 {" Ka1 = [1 0];
6 k% I) ?) F( x. ]- _" g3 y; \dev1 = [0.01 10^-8];( H# g/ [' q# n. z' y
fs1 = fs*L1;
9 p3 p; S' ~/ _8 |+ g4 r[n1 fo1 ao1 w1] = firpmord(f1,a1,dev1,fs1);
) h( f4 A. K% b- fb1 = firpm(n1,fo1,ao1,w1);

%%Second Filter3 b) J+ `5 c  J, J' ]
f2 = [20160 30240];
& c4 C; c6 g1 p  J  ra2 = [1 0];
8 _) O, y. O1 Adev2 = [0.01 10^-8];
! B+ W$ d( E; t. ~. bfs2 = fs1*L2/M1;
, B& t/ @- F, S) p$ R( H[n2 fo2 ao2 w2] = firpmord(f2,a2,dev2,fs2);
7 f6 R4 k3 B6 _9 Vb2 = firpm(n2,fo2,ao2,w2);

%%Third Filter
% l- O3 g! V7 f( S$ _8 O6 Ff3 = [16800 50400];0 J7 r& Z6 [. U  @# e
a3 = [1 0];
% l1 ]$ O. o3 p, W9 _' v5 B! _dev3 = [0.01 10^-8];7 b, j, J( W4 \5 {
fs3 = fs2*L3/M2;
; y  T7 u3 V3 t8 h" D# p$ \+ \' T$ J[n3 fo3 ao3 w3] = firpmord(f3,a3,dev3,fs3);: \# H" _. n8 u
b3 = firpm(n3,fo3,ao3,w3);