mulitplex1.m

多通道滤波器组的程序

tic; %beginning time
clear; clc;
%%
%====================
M = 3; % there are 3 channels
N = 256;
n = -N : N;
fs = 100;  % sample frequency
dt = 1/fs;
t = n * dt;
t1 = linspace(-N*dt,N*dt,((2*N+1)*M));
mu0 = -20;   % the chirp frequency
f0 = 5;   % center frequency of the first signal x0
f1 = 35;   % center frequency of the second signal x1
f2 = 70;   % center frequency of the third signal x2
alfa = acot(-mu0);
du=2*pi*abs(sin(alfa))/(dt*(2*N+1));
u = linspace(-pi*abs(sin(alfa)) , pi*abs(sin(alfa)) , 2 * N + 1);
u1 = linspace(-pi*abs(sin(alfa)) , pi*abs(sin(alfa)) ,((2*N+1)*M));
uu = linspace(0 , 2*pi*abs(sin(alfa)) , 2 * N + 1);
uu1 = linspace(0 , 2*pi*abs(sin(alfa)) ,((2*N+1)*M));

%%
x0 = exp(j * 0.5 * mu0 * t .^ 2  + j * 2 * pi * f0 * t);
x1 = exp(j * 0.5 * mu0 * t .^ 2  + j * 2 * pi * f1 * t);
x2 = exp(j * 0.5 * mu0 * t .^ 2  + j * 2 * pi * f2 * t);

%%%%% the frft of the original signal 
xx0 = dfrft3(x0 , alfa , N , N , dt);
xx1 = dfrft3(x1 , alfa , N , N , dt);
xx2 = dfrft3(x2 , alfa , N , N , dt);
%%
max0_1 = max(abs(xx0));
max1_1 = max(abs(xx1));
max2_1 = max(abs(xx2));
%%
%%%%% to plot the original signal and its frft sprctrum
%%%========= for the first signal
figure(1); 
subplot(211);
plot(t,real(x0));
title('the first original signal x0');
ylabel('x0(t)');
xlabel('t')
% axis([-N*dt N*dt -1 1]);

subplot(212);
plot(u,abs(fftshift(xx0)));
title('the frft of the first original signal x0');
ylabel('xx0(u)');
xlabel('u');
%%%========= for the second signal
figure(11); 
subplot(211);
plot(t,real(x1));
title('the second original signal x1');
ylabel('x1(t)');
xlabel('t')
% axis([-N*dt N*dt -1 1]);

subplot(212);
plot(u,abs(fftshift(xx1)));
title('the frft of the second original signal x1');
ylabel('xx1(u)');
xlabel('u');
%%%========= for the third signal
figure(12); 
subplot(211);
plot(t,real(x2));
title('the third original signal x2');
ylabel('x1(t)');
xlabel('t')
% axis([-N*dt N*dt -1 1]);

subplot(212);
plot(u,abs(fftshift(xx2)));
title('the frft of the third original signal x2');
ylabel('xx1(u)');
xlabel('u');
%%
%%%%%% interpolation by M = 3
%%%
N2 = floor(M*(2*N+1)-1)/2;

%%%s0 is the first interpolated signal
s0 = upsample( x0 , M);
ss0 = dfrft3(s0 , alfa , N2 , N2 , dt/M); %ss0 is the first frft of s0 with angle of alfa
max0_2 = max(abs(ss0));

%%%s1 is the second interpolated signal
s1 = upsample( x1 , M);
ss1 = dfrft3(s1 , alfa , N2 , N2 , dt/M); %ss1 is the second frft of s0 with angle of alfa
max1_2 = max(abs(ss1));

%%%s2 is the third interpolated signal
s2 = upsample( x2 , M);
ss2 = dfrft3(s2 , alfa , N2 , N2 , dt/M); %ss2 is the third frft of s0 with angle of alfa
max2_2 = max(abs(ss2));

%%
%%%%%% filtering
h0 = zeros(1 , 2 * N * M + 3);
h1 = h0;
h2 = h0;

%%%======== filtering for the first interpolated signal
h0([2*N+2 : 4*N+2]) = M*ones(1, length([2*N+1 : 4*N+1])); % the first filter 
% ss0 = fftshift(ss0);
v0 = ss0 .* h0;
% v0 = ifftshift(v0);

%%%======== filtering for the second interpolated signal
h1([4*N+3 : 6*N+3]) = M*ones(1, length([4*N+2 : 6*N+2])); % the second filter 
% ss1 = fftshift(ss1);
v1 = ss1 .* h1;
% v1 = ifftshift(v1);

%%%======== filtering for the third interpolated signal
h2([1 : 2*N+1]) = M*ones(1, length([1 : 2*N+1])); % the third filter 
% ss2 = fftshift(ss2);
v2 = ss2 .* h2;
% v2 = ifftshift(v2);

%%
%%%%%% to plot the spectrum which is the frft of the signal after
%%%%%% interpolating the original signal
%%%======== for the first signal
figure(2);
subplot(211);
plot(u1,abs((ss0))/max0_2 * max0_1);
title('the frft of the first interpoltated signal for x0');
ylabel('ss0(u)');
xlabel('u');

subplot(212);
plot(u1,abs((v0))/max0_2 * max0_1);
title('the spectrum v0 which is the result of filtering ss0 ');
ylabel('v0(u)');
xlabel('u');

%%%======== for the second signal
figure(21);
subplot(211);
plot(u1,abs((ss1))/max1_2 * max1_1);
title('the frft of the second interpoltated signal for x1');
ylabel('ss1(u)');
xlabel('u');

subplot(212);
plot(u1,abs((v1))/max1_2 * max1_1);
title('the spectrum v1 which is the result of filtering ss1 ');
ylabel('v1(u)');
xlabel('u');

%%%======== for the third signal
figure(22);
subplot(211);
plot(u1,abs((ss2))/max2_2 * max2_1);
title('the frft of the third interpoltated signal for x2');
ylabel('ss2(u)');
xlabel('u');

subplot(212);
plot(u1,abs((v2))/max2_2 * max2_1);
title('the spectrum v1 which is the result of filtering ss2 ');
ylabel('v2(u)');
xlabel('u');

%%
%%%%%% inverse frft of the spectrum after being filtered to get the
%%%%%% interpolated signal in the time domain

iv0 = dfrft3(v0 , -alfa , N2 , N2 ,  du);
iv1 = dfrft3(v1 , -alfa , N2 , N2 , du);
iv2 = dfrft3(v2 , -alfa , N2 , N2 , du);

%%
%%%======== to plot the signal which is mixed with the three interpplated
%%%signals s0 , s1, s2

figure(5);
subplot(211);
plot(u1 ,abs((v0)) , u1 , abs((v1)) , 'r' , u1 , abs((v2)) , 'g');
legend ('the first spectrum' , 'the second spectrum' , 'the third spectrum');
title ('the mixed spectrum v(u)');
ylabel('v(u)');
xlabel('u');

subplot(212);
plot(t1 , real(iv0) + real(iv1) + real(iv2));
title('the mixed signal v(t)');
ylabel('v(t)');
xlabel('t');

%%
%%%%% decimation form the interpolated signals in the time doamin
y0 = downsample(iv0 , M);
y1 = downsample(iv1 , M);
y2 = downsample(iv2 , M);

%%
%%%======== for the first signal x0
figure(3);
subplot(211);
plot(t1,real(iv0));
title('the first interpolated signal after being filtered for x0');
ylabel('iv0(n)');
xlabel('t');
%  axis([-N*dt N*dt -1 1]);

subplot(212);
plot(t,real(y0));
title('the first recovered signal of x0');
ylabel('y0');
xlabel('t');
% axis([-N*dt N*dt -1 1]);

%%%======== for the second signal x1
figure(31);
subplot(211);
plot(t1 , real(iv1));
title('the second interpolated signal after being filtered for x1');
ylabel('iv1(n)');
xlabel('t');

subplot(212);
plot(t , real(y1));
title('the second recovered signal of x1');
ylabel('y1');
xlabel('t');

%%%======== for the third signal x2
figure(32);
subplot(211);
plot(t1 , real(iv2));
title('the third interpolated signal after being filtered for x2');
ylabel('iv2(n)');
xlabel('t');

subplot(212);
plot(t , real(y2));
title('the third recovered signal of x2');
ylabel('y2');
xlabel('t');

%%
%%%%% to compare the signals between the originals and the recovered ones
%%%======== for the first signal
figure(4);
subplot(211);
plot(t,real(y0),t, real(x0) , 'r:');
title('the first original signal x0 and the first recovered signal y0');
ylabel('y0 and x0');
xlabel('t');
legend('the recovered signal y0' ,'the original signal x0');
% axis([-N*dt N*dt -1 1]);

subplot(212);
plot(real(y0) - real(x0));
title('errors between the first signals x0 and y0');
axis([0 520 -10^(-13) 10^(-13)]);

%%%======== for the second signal
figure(41);
subplot(211);
plot(t , real(y1) , t ,real(x1), 'r:');
title('the second original signal x0 and the first recovered signal y1');
ylabel('y1 and x1');
xlabel('t');
legend('the recovered signal y1' ,'the original signal x1');

subplot(212);
plot(real(y1) - real(x1));
title('errors between the second signals x1 and y1');
axis([0 520 -10^(-13) 10^(-13)]);

%%%======== for the third signal x2
figure(42);
subplot(211);
plot(t , real(y2) , t ,real(x2), 'r:');
title('the third original signal x0 and the first recovered signal y2');
ylabel('y2 and x2');
xlabel('t');
legend('the recovered signal y2' ,'the original signal x2');

subplot(212);
plot(real(y2) - real(x2));
title('errors between the third signals x2 and y2');
axis([0 520 -10^(-13) 10^(-13)]);

%======================
toc; % end time