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

loveeatmore

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

%input data
fa = 8000; %%signal frequency
% s0 P' R( @' }+ l3 V' D8 `fs = 44100; %%44.1kHz sampling frequency
# Q; R% T) P; j- bn = 1:64;
  c* ^9 T5 @" W( P) Q% v& Wx = sin(2*pi*fa*n/fs);
lengthx = length(x);
+ v; d8 s) J5 H. A' O. E9 z, |" Et = n*1000/128/fs;
%plot(t,x);xlabel('time in ms');

% [y fs nbits] = wavread('acoustic.wav');
" W, Y7 w  k2 Z3 [% x = y';
% lx = length(x);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ n9 [+ d1 u. P' ^) r! o( D, _7 z%Resample parameters
2 B9 p& U* C1 v# S; W. y6 c( CL1 = 8;
: {* d; y4 }: X& C4 k6 }/ SM1 = 7;
7 w4 v! G( K" L+ }8 A/ hL2 = 4;
M2 = 3;
" i3 k- S2 b! d- j3 I" RL3 = 10;
, D7 \0 `3 d/ d+ ^M3 = 7;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%filters
load filter.mat
filter1 = b1;
6 c6 b6 U* i1 }- U1 tfilter2 = b2;
, [) G6 D+ N% n+ b9 r4 Lfilter3 = b3;
9 Y! f6 K" ~& p3 z+ `2 d  A%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Initialize the subfilters

%First Filter
% G# b; g' h/ D# \2 w8 d/ TupFilter1 = decompose(filter1,L1);
8 v3 L/ Y' P3 L7 ^/ X" i5 udownFilter1 = decompose(filter1,M1);

%Second Filter
upFilter2 = decompose(filter2,L2);
% G8 ^5 c( r9 q' gdownFilter2 = decompose(filter2,M2);

%Second Filter
8 L7 z5 A+ E8 h/ g2 x3 UupFilter3 = decompose(filter3,L3);
downFilter3 = decompose(filter3,M3);

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

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
6 y) Y) T  ^3 a4 q3 M, k1 f%Polyphase Filtering

%First Filter
3 i( p) r- h$ q7 U%
%%upsampling
intermediatex1 = zeros(L1,length(x));
) X& A, k0 ^" H/ vx11 = zeros(L1,length(x)*L1);
, Y9 Z6 L- o5 G+ I  M$ `2 A  r  gfor i = 11
intermediatex1(i, = filter(upFilter1(i,,1,x); %%do the filtering before upsampling
x11(i, = upsample(intermediatex1(i,:),L1); %%upsample by L, inserting L-1 zeros.
$ j9 x8 s( ]* p) ^1 o6 ^end

X1 = zeros(1,length(x)*L1);
for i = 1:length(x)
X1(L1*(i-1)+11*(i-1)+L1) = x11(:,(i-1)*L1+1)';
9 B, u5 f  L" ?* b9 vend
%%downsampling
X1_down = [X1 zeros(1,M1-mod(length(X1),M1))];
output1 = zeros(1,length(X1_down)/M1);
x1down = zeros(M1,length(X1_down)/M1);
x11down = x1down;
for i = 1:M1
0 {5 C( J- m0 v1 T" rfor k = 1:length(X1_down)/M1
3 H; }# ]0 \2 i8 L+ Bx1down(i,k) = X1_down(M1*(k-1)+i);
end
& B6 {5 s, _1 Kx11down(i,:) = filter(downFilter1(1,:),1,x1down(i,:));
output1 = output1 + x11down(i,:);
5 y2 z# d, w6 q, j+ lend

%Second Filter
- Y' H( Y9 i6 f$ C; x. G%
%%upsampling
3 T; V: r  P) p4 P! vintermediatex2 = zeros(L2,length(output1));
x22 = zeros(L2,length(output1)*L2);
! D! x2 o' B" j: F6 o: F4 v5 z' ]for i = 12
intermediatex2(i,:) = filter(upFilter2(i,:),1,output1); %%do the filtering before upsampling
x22(i,:) = upsample(intermediatex2(i,:),L2); %%upsample by L, inserting L-1 zeros.
, O8 m& h& Y7 H2 Jend

X2 = zeros(1,length(output1)*L2);
: r+ t  \7 i* k7 @  `! l4 F, Sfor i = 1:length(output1)
  I& u( ~, K* w8 z, a1 ?  [/ R3 {X2(L2*(i-1)+1:L2*(i-1)+L2) = x22(:,(i-1)*L2+1)';
# k9 n2 L% @" T7 V/ _# iend
%%downsampling
X2_down = [X2 zeros(1,M2-mod(length(X2),M2))];
output2 = zeros(1,length(X2_down)/M2);
8 `6 e7 I( U# h2 y6 Xx2down = zeros(M2,length(X2_down)/M2);
x22down = x2down;
* m$ o+ Z2 o8 i8 g! Xfor i = 1:M2
. U  ^$ \6 E( h3 k% o. ~& O/ @for k = 1:length(X2_down)/M2
x2down(i,k) = X2_down(M2*(k-1)+i);
end
1 p* @. r1 F, b* X0 f8 ~5 |x22down(i,:) = filter(downFilter2(1,:),1,x2down(i,:));
output2 = output2 + x22down(i,:);
6 H* c8 q3 O8 L0 K0 c* `end

%Third Filter
%
7 H' c; o$ O. U4 \. j# _- n1 {8 u%%upsampling
# I7 V% E: ^5 uintermediatex3 = zeros(L3,length(output2));
9 g& s6 d" M4 S$ A, z/ Yx33 = zeros(L3,length(output2)*L3);
for i = 1:L3
intermediatex3(i,:) = filter(upFilter3(i,:),1,output2); %%do the filtering before upsampling
# g9 s+ P  r, wx33(i,:) = upsample(intermediatex3(i,:),L3); %%upsample by L, inserting L-1 zeros.
end

X3 = zeros(1,length(output2)*L3);
" `0 Y! o; v  Dfor i = 1:length(output2)
# K: a( Y& a( ~X3(L3*(i-1)+1:L3*(i-1)+L3) = x33(:,(i-1)*L3+1)';
1 ^; d2 }% h3 k6 Xend
8 E( r8 E- Z$ b" Q9 X# r%%downsampling
$ ~+ @' Y2 C3 E* \X3_down = [X3 zeros(1,M3-mod(length(X3),M3))];
output3 = zeros(1,length(X3_down)/M3);
x3down = zeros(M3,length(X3_down)/M3);
x33down = x3down;
for i = 1:M3
3 Y* F" n( l$ a$ d+ t' {, sfor k = 1:length(X3_down)/M3
) A4 |0 K  I; h- a/ W, cx3down(i,k) = X3_down(M3*(k-1)+i);
% p& U+ u5 h/ i2 B4 rend
x33down(i,:) = filter(downFilter3(1,:),1,x3down(i,:));
6 m0 t& J0 c/ U- F+ u* P1 youtput3 = output3 + x33down(i,:);
end

output3 = output3/max(output3);
lout = length(output3);

nx = 1:length(output3);
7 F1 o  w  f3 p. C/ l6 |& Rxorigin = sin(2*pi*fa*nx/96000);
error1 = xorigin - output3;
yin = [0 0 resample(x,320,147)];
) L7 t1 Q1 `/ q9 _/ t. n- ]0 Jerror2 = xorigin - yin;

figure(1)
7 \! y" b/ M5 }3 u9 P. A  o0 }subplot(3,1,1);plot(xorigin);title('True Signal');
subplot(3,1,2);plot(output3);title('Resampled Signal');
subplot(3,1,3);plot(error1);title('Error');

figure(2)
% F* J- d. e% }( psubplot(3,1,1);plot(xorigin);title('True Signal');
subplot(3,1,2);plot(yin);title('Built-in Resampled Signal');
/ x- B$ S. K1 v. \/ K% \2 I4 b: Csubplot(3,1,3);plot(error2);title('Error');

---

function x_decompose = decompose(x,factor)
lx = length(x);

if mod(lx,factor) ~= 0;
9 V  s6 m$ w1 b7 Y7 RX = [x zeros(1,factor-mod(lx,factor))];
0 C$ H8 U' j2 G. P8 J, X) e; aelse
X = x;
end

x_decompose = zeros(factor,length(X)/factor);
for i = 1:factor
, @6 c# }  l- u" |for k = 1:length(X)/factor
$ l+ x1 P6 @; y  X* ex_decompose(i,k) = X(factor*(k-1)+i);
end
0 n5 w; g# Z* X* D. e6 f# L2 Bend

end

%Filter Generation
fs = 44100;
2 M( ~& i8 j! c: lL1 = 8;
2 K2 x- R* w0 L, y2 Z$ y% OM1 = 7;
L2 = 4;
4 r" r) A, Z: @! A. G. C, S' O" ^* sM2 = 3;
4 F6 {" a0 y1 ^3 K$ K( N  {& C+ l6 W& {L3 = 10;
) r9 G' K5 e! N1 sM3 = 7;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
- Q( ~0 g; o( X! _! t, I%
%%First Filter
f1 = [20000 24000];
a1 = [1 0];
8 [! B: T0 h) c& D" q. Y/ T5 l, `dev1 = [0.01 10^-8];
" T3 Z9 v! o0 |$ J* f( J0 vfs1 = fs*L1;
3 b6 W+ {3 `& i/ A- j1 t[n1 fo1 ao1 w1] = firpmord(f1,a1,dev1,fs1);
& g. p2 |# g. W1 v9 o  Bb1 = firpm(n1,fo1,ao1,w1);

%%Second Filter
f2 = [20160 30240];
9 B4 @! u9 i# N8 J8 Ba2 = [1 0];
dev2 = [0.01 10^-8];
+ L$ I* w) r8 y& l$ b& qfs2 = fs1*L2/M1;
( o& W. h& f% u* @' y5 o$ E[n2 fo2 ao2 w2] = firpmord(f2,a2,dev2,fs2);
9 U: s. g6 ^; W; Q1 C0 Ub2 = firpm(n2,fo2,ao2,w2);

%%Third Filter
- Z+ B5 s+ H( G' ?' ^, zf3 = [16800 50400];
a3 = [1 0];
5 N+ U6 g; Q% G; K; V; Bdev3 = [0.01 10^-8];
fs3 = fs2*L3/M2;
[n3 fo3 ao3 w3] = firpmord(f3,a3,dev3,fs3);
b3 = firpm(n3,fo3,ao3,w3);

木棉花_MM

这个……