以下是使用相关 32 通道分析和合成滤波器组的 32 分量多相分解的工作代码。正如我已经评论过的，由于信号和滤波器长度较短，这种 CAE 的速度增益并不显着。然而，进一步的架构改进以及编码优化可以提供更好的结果。

% S0 - Load the prototype lowpass filter impulse response h0[n]:
% --------------------------------------------------------------
load h2.mat;         % h[n] is the prototype lowpass filter of length 512
L = length(h);

% S1 - Create the 32 x 512 analysis filter-bank hha[k,n] by cosine modulation from protoype :
% -----------------------------------------------------------------------------------------
numbands = 32;                 % number of banks (channels)
n=0:L-1;

hha=zeros(numbands,L);         % bank of filters hha[k,n] = 32 x 512 array.
for k=0:1:numbands-1 
   hha(k+1,:) = h.*cos( ( (2*k+1)*pi*(n-16) ) / (2*numbands) );
end


% S2 - Create the 32-polyphase components hhap[k,m,n] , for each one of 32 analysis filters hha[k,n]:
% ---------------------------------------------------------------------------------------------------
numpoly = numbands;             % polyphase component number = decimation ratio = number of channels
hhap = zeros(numbands,numpoly, L/numpoly);  % hhap = 32 x 32 x 512/32 , 3D ANALYSIS filter bank array

M = numpoly;                    % polyphase system decimation ratio
for k=1:numbands
    for m = 1:numpoly
        hhap(k,m,:) = hha(k,m:M:end);       % create the m-th polyphase component of k-th channel filter
    end
end


% S3 - Design the 32 x 512  synthesis (complementary) filter bank :
% -----------------------------------------------------------------
numbands = 32;                  % number of banks
n=0:L-1;
hhs = zeros(numbands,L);        % bank of filters
for k=0:1:numbands-1 
   hhs(k+1,:) = h.*cos( ( (2*k+1)*pi*(n+16) ) / (2*numbands) );
end


% S4 - Obtain the 32-polyphase components hhsp[k,m,n] , for each one of 32 synthesis filters hhs[k,n]:
% ----------------------------------------------------------------------------------------------------
numpoly = numbands;             % polyphase component number = interpolation ratio = number of channels
hhsp = zeros(numbands,numpoly, L/numpoly);  % hhap = 32 x 32 x 512/32 , 3D ANALYSIS filter bank array
M = numpoly;                    % polyphase system decimation ratio
for k=1:numbands
    for m = 1:numpoly
        hhsp(k,m,:) = hhs(k,m:M:end);       % create the m-th polyphase component of k-th channel filter
    end
end


% S5 - Generate the test input signal
% -----------------------------------
N = 2*1024;
wav_in = cos(0.01791*pi*[0:N-1]);        % pure sine tone

% S6 - Apply test signal to the filterbank: ANALYSIS STAGE :
% -----------------------------------------------------------
yyd = zeros( numbands, floor(N/numbands));   % decimated outputs..
M = numbands;
for k=1:1:numbands       
    temp = conv([wav_in(1:M:end),0] , hhap(k,1,:));
    for m=2:M
        temp = temp + conv([0,wav_in(M-m+2:M:end)],hhap(k,m,:));   
    end
    yyd(k,:) = temp(L/(2*M)+1 : L/(2*M)+N/numbands);
end

% S7 - Apply SYNTHESIS filterbank on the decimated signal :
% ---------------------------------------------------------
ys = zeros(1, N);

for k=1:numbands
    temp = zeros(1, N+L-1);
    for m = 1:numpoly
        temp(m:numbands:end-31) = conv( yyd(k,:) , hhsp(k,m,:) );
    end

    ys = ys + temp(L/2+1:L/2+N);    
end
ys = numbands*ys;


% SX - DISPLAY RESULTS:
% ---------------------
L = length(h);
figure,subplot(2,1,1)
stem([0:L-1],h);title('The Prototype Lowpass Filter');
subplot(2,1,2)
plot(linspace(-1,1,4*L),20*log10(abs(fftshift(fft(h,4*L)))));
grid on;

figure
plot(linspace(-1,1,4*L),20*log10(abs(fftshift(fft(hha(1,:),4*L)))));
hold on
for k=2:numbands
    plot(linspace(-1,1,4*L),20*log10(abs(fftshift(fft(hha(k,:),4*L)))));
end
title('32 CHANNEL FILTERBANK');

figure,subplot(2,1,1)
plot(wav_in);title('input signal')
subplot(2,1,2)
plot(linspace(-1,1,4*N),20*log10(abs(fftshift(fft(wav_in,4*N)))));

figure,subplot(2,1,1)
plot(ys);title('Synthesized Back');
subplot(2,1,2)
plot(linspace(-1,1,4*N),20*log10(abs(fftshift(fft(ys,4*N)))));