В результаті обчислень отримуємо передаточну функцію смугового фільтра Чебишева І роду:
Рис. 7. діаграма нулів та полюсів передаточної функції фільтра.
Рис. 5. АЧХ смугового фільтра Чебишева І роду.
Рис. 6. ФЧХ смугового фільтра Чебишева І роду.
Рис. 7. АЧХ фільтра в порівнянні з бажаною АЧХ.
Щоб перевірити правильність розрахунків промоделюємо роботу фільтра у середовищі MatLab. Для моделювання використаємо тестовий сигнал що складається з трьох гармонічних сигналів різної частоти, фази та амплітуди.
Рис. 8. Вхідний та вихідний сигнали.
1.5 Побудова спектрів вхідного та вихідного сигналів
Побудуємо спектри вхідного та вихідного сигналів. Для побудови спектрів використаємо алгоритм неперервного перетворення Фур’є.
Для дійсних функцій ряд Фур’є часто означається, як ряд по синусах та косинусах, аргументи яких обираються в такий спосіб, аби на відрізку означення функції вкладалась ціла кількість періодів. Легко зрозуміти, що перелік таких синусів та косинусів вичерпується множинами
та 
, де n є довільним натуральним числом, а для 
може приймати також значення 0. 
%------------------------- Побудова спектрів ---------------------------%
[y,x]=lsim(Hrez,Analog_sign,t);
a0=dt*sum(Analog_sign)/T1;
for i=1:k
ak(i)=(2*dt/T1)*sum(Analog_sign.*cos(t*i*W1));
bk(i)=(2*dt/T1)*sum(Analog_sign.*sin(t*i*W1));
end;
Ak=sqrt(ak.^2+bk.^2);
k1=1:i; % номера відліків
W0=k1*W1; % масштабовані відліки
Wn=zeros(1,50); % завдання розміру вектора відліків
Wn=0:W0:50*W0; % заповнення вектора відліків (перші 50 елементів вектора W0)
Akn=zeros(1,50); % завдання розміру вектора амплітуд з урахуванням а0
Akn(1,1:50)=Ak(1,1:50); %заполнение вектора амплитуд
out_sign=y';
for i=1:k
ak_out(i)=(2*dt/T1)*sum(out_sign.*cos(t*i*W1));
bk_out(i)=(2*dt/T1)*sum(out_sign.*sin(t*i*W1));
end;
a0_out=dt*sum(out_sign)/T1;
Ak_out=sqrt(ak_out.^2+bk_out.^2);
Akn_out=zeros(1,50);% завдання розміру вектора амплітуд з урахуванням a0
Akn_out(1,1:50)=Ak_out(1,1:50);%заполнение вектора амплитуд
figure;
hold on
stem(Wn,[a0 Akn],'ro');
stem(Wn,[a0_out Akn_out],'b+');
xlabel('w'); ylabel('Magnitude');
axis([0 150000 0 70]);
hold off
grid on
Рис. 9. Спектри вхідного та вихідного сигналу.
2. Розрахунок та дослідження цифрового фільтра
2.1 Визначення порядку та частоти зрізу цифрового фільтра
Як і у випадку проектування аналогових фільтрів для визначення порядку цифрових фільтрів та їх частоти зрізу можна скористатися функціями пакета Communication середовища MatLab.
Згідно з теоремою Котельникова частота дискретизації має бути принаймні вдвічі більшою ніж максимальна частота сигналу. Для моделі фільтра використаємо частоту
, де 
- частота дискретизації, 
- верхня межи смуги затримки.% Початкові дані
Wp1=8442; % Нижня межа смуги пропускання рад / сек
Wp2=12940; % Верхня межа смуги пропускання рад / сек
Ws1=7914; % Нижня межа смуги затримки рад / сек
Ws2=13468; % Верхня межа смуги затримки рад / сек
Rp=2; % Коефіцієнт пульсацій дБ
Rs=37; %Затухання дБ
Wd=Ws2*4 % Частота дискретизації рад / сек
Wpz=[Wp1/Wd Wp2/Wd]
Wsz=[Ws1/Wd Ws2/Wd]
[n,wn]=cheb1ord (Wpz, Wsz, Rp, Rs);
В результаті розрахунків маємо:
Мінімальний порядок фільтра n=9.
Частоти зрізу (відносні) Wn = [0.1567 0.2402]
2.2 Визначення коефіцієнтів фільтра та отримання передаточної функції
Для визначення коефіцієнтів фільтра використаємо середовище MatLab, при цьому можливості MatLab дозволяють виключити етап розрахунку фільтра прототипу і одразу визначити коефіцієнти фільтра згідно з ТЗ.
Wp1z=Wp1/Wd
Wp2z=Wp2/Wd
[b,a]=cheby1(n,Rp,[Wp1z Wp2z]);
Hrz=tf(b,a,1/(Wd*2*pi))
В результаті розрахунків отримуємо:
Вектор коефіцієнтів чисельника
b = 1.0e-008 * [0.0053 0 -0.0480 0 0.1920 0 -0.4481 0 0.6721 0 -0.6721 0 0.4481 0 -0.1920 0 0.0480 0 -0.0053]
a = 1.0e+004 * [0.0001 -0.0014 0.0102 -0.0462 0.1510 -0.3780 0.7496 -1.2033 1.5843 -1.7234 1.5525 -1.1555 0.7054 -0.3486 0.1365 -0.0409 0.0088 -0.0012 0.0001]
Передаточна функція фільтра:
Побудуємо характеристики отриманого фільтра за допомогою MatLab.
[b,a]=cheby1(n,Rp,[Wp1z Wp2z]);
[H,W]=freqz(b,a,'whole');
figure;
plot (W,abs(H));grid on;
Hrz=tf(b,a,1/(Wd*2*pi))
figure;
pzmap(Hrz)
Рис. 10 Діаграма нулів та полюсів фільтра.
Рис. 11 АЧХ смугового фільтра.
Рис. 12. ФЧХ смугового фільтра.
2.3 Моделювання роботи дискретного фільтра
Щоб перевірити правильність розрахунків промоделюємо роботу дискретного фільтра у середовищі MatLab. Для моделювання використаємо тестовий сигнал що складається з трьох гармонічних сигналів різної частоти, фази та амплітуди.
%-----------------Формирование исходного сигнала--------------------------%
k1=10;
k2=9;
k3=12;
phi1=pi;
phi2=pi/4;
phi3=pi/2;
W1=6000;
W2=10976;
W3=24387;
sign_freq=[W1 W2 W3];
sign_phase=[phi1 phi2 phi3];
sign_koeff=[k1 k2 k3];
Fs=Wd/(2*pi);
td=0:1/Fs:0.1;
Signal_digital=0;
for i=1:3
Signal_digital=Signal_digital+sign_koeff(i).*sin(sign_freq(i).*td);
end;
wdd=2*pi*Fs;
wd_mod=0:wdd/2/15671:wdd/2;
Sd_in=fft(Signal_digital);
Signd_out=filter(b,a,Signal_digital);
Sd_out=fft(Signd_out);
figure
hold on
plot(td(1,1:700), Signal_digital(1,1:700),'g-');
plot(td(1,1:700), Signd_out(1,1:700),'b-')
hold off;
grid on;
axis([0 0.08 -35 35]);
xlabel('t, c');
ylabel('Amplitude');
title ('Filtered and NOT filtered signals');
Рис. 13. Вхідний та вихідний сигнали.
2.4 Побудова спектрів вхідного та вихідного сигналів
Побудуємо спектри вхідного та вихідного сигналів. Для побудови ряду використаємо дискретне перетворення Фур’є. Припустимо, що у нас задано функцію дискретної змінної на рівномірній гратці, тобто задані значення функції fk для скінченної послідовності значень аргументу xk = kh – таблиця функції {fk; kh; k = 0::N} . Тут за допомогою
позначено крок гратки – відстань між сусідніми вузлами. Перетворення Фур’є такої функції можна означити, як суму
