clear all;
close all;
clc;
clf;

% Kazda prenosova funkce ma svuj vektor frekvence kvùli nezavislosti
p =  j*(10:10000:100000000); 
p1 = j*(10:10000:100000000); 
p2 = j*(10:10000:100000000); 
p3 = j*(10:10000:100000000); 
p4 = j*(10:10000:100000000); 
p5 = j*(10:10000:100000000); 

% A) Tolerantcni schema zprisnene co se tyka frekvence amin 0.2 amax 51

%---------------------------------- Prenosove funkce filtru -------------------------------------
%----------------------------------- Cauerova aproximace B6 -------------------------------------
H = (0.2989442e+14.*(p.^2 + 0.2641934e+16).*(p.^2 + 0.1557567e+16))./...
((p.^2 + 0.1327380e+8.*p + 0.6081370e+15).*(p.^2 + 3657890.*p + 0.8485143e+15).*...
(p.^2 + 0.2423799e+8.*p + 0.2439486e+15));
%----------------------------------- Cauerova aproximace A6 -------------------------------------
H1 = (0.2609761e-2.*(p1.^2 + 0.2346993e+16).*(p1.^2 + 0.1481169e+16).*(p1.^2 + 0.1463616e+17))./...
((p1.^2 + 0.1287481e+8.*p1 + 0.6199607e+15).*(p1.^2 + 3442557.*p1 + 0.8462264e+15).*...
(p1.^2 + 0.2479755e+8.*p1 + 0.2589963e+15));
%---------------------------------- Chebisevova aproximace 11 -----------------------------------
H2= .4245405e80./...
((p2 + 5788554).*(p2.^2 + 0.1110815e+8.*p2 + 0.9721541e+14).*(p2.^2 + 9739284.*p2 + 0.2681127e+15).*...
(p2.^2 + 7581394.*p2 + 0.4919405e+15).*(p2.^2 + 4809305.*p2 + 0.6976352e+15).*...
(p2.^2 + 1647594.*p2 + 0.8198900e+15));
%--------------------------------- Butterworthova aproximace 26 ---------------------------------
% Nelze zobrazit skupinové zpozdeni
 H3 = 0.7488257e+202./...
 ((p3 + 0.2997981e+8).*(p3.^2 + 0.5955420e+8.*p3 + 0.8987892e+15).*...
 (p3.^2 + 0.5834341e+8.*p3 + 0.8987892e+15).*...
 (p3.^2 + 0.5634362e+8.*p3 + 0.8987892e+15).*(p3.^2 + 0.5358188e+8.*p3 + 0.8987892e+15).*...
 (p3.^2 + 0.5009554e+8.*p3 + 0.8987892e+15).*(p3.^2 + 0.4593174e+8.*p3 + 0.8987892e+15).*...
 (p3.^2 + 0.4114679e+8.*p3 + 0.8987892e+15).*(p3.^2 + 0.3580541e+8.*p3 + 0.8987892e+15).*...
 (p3.^2 + 0.2997981e+8.*p3 + 0.8987892e+15).*(p3.^2 + 0.2374879e+8.*p3 + 0.8987892e+15).*...
 (p3.^2 + 0.1719661e+8.*p3 + 0.8987892e+15).*(p3.^2 + 0.1041188e+8.*p3 + 0.8987892e+15).*...
 (p3.^2 + 3486342.*p3 + 0.8987892e+15));
%----------------------------- Inverzní Chebisevova aproximace A11 ------------------------------
H4= (646478.8.*(p4.^2 + 0.1783403e+17).*(p4.^2 + 0.4842904e+16).*(p4.^2 + 0.2478379e+16).*...%
		(p4.^2 + 0.1710772e+16).*(p4.^2 + 0.1444807e+16))./...
      ((p4 + 0.5395709e+8).*(p4.^2 + 0.8901189e+8.*p4 + 0.2502792e+16).*...
      (p4.^2 + 0.5669833e+8.*p4 + 0.1818285e+16).*(p4.^2 + 0.3249577e+8.*p4 + 0.1338741e+16).*...
      (p4.^2 + 0.1659242e+8.*p4 + 0.1077571e+16).*(p4.^2 + 5093699.*p4 + 0.9656097e+15));
%----------------------------- Inverzní Chebisevova aproximace B12 ------------------------------
H5= (0.4611081e+14.*(p5.^2 + 0.1083282e+17).*(p5.^2 + 0.3965087e+16).*(p5.^2 + 0.2289244e+16).*...
   (p5.^2 + 0.1675843e+16).*(p5.^2 + 0.1451323e+16))./...
   ((p5.^2 + 0.7574424e+8.*p5 + .2129478e+16).*(p5.^2 + .4893114e+8.*p5 + 0.1607206e+16).*...
    (p5.^2 + 0.2906651e+8.*p5 + .1246967e+16).*(p5.^2 + .1525892e+8.*p5 + 0.1042633e+16).*...
    (p5.^2 + 4750408.*p5 + 0.9521779e+15).*(p5.^2 + 0.9965325e+8.*p5 + 0.2602779e+16));
  


a = 20*log10(abs(H));		% Cauer B6
a1 = 20*log10(abs(H1));		% Cauer A6  *
a2 = 20*log10(abs(H2));		% Chebysev 11  *
a3 = 20*log10(abs(H3));		% Butter 26   *
a4 = 20*log10(abs(H4));		% ICheby A11 *
a5 = 20*log10(abs(H5));		% ICheby B12 


%---------------- Porovnání vysledných prenosovych funkci jednotlivych aproximaci ---------------


% Zobrazeni propustneho pasma 
figure(1)
semilogx(imag(p)/(2*pi),a1,'k');
hold on;
text(3500000,-0.1,'CauerA6'); hold on;
axis([3000000 5000000 -2 1]);
xlabel('--> f[Hz]');
ylabel('--> a[dB]');
title('Zobrazení propustného pásma filtru');
hold on;
% Zobrazeni tolerancniho schematu
semilogx([4500000 4500000],[-80 -0.2],'k'); hold on;
semilogx([6000000 6000000],[-51 0],'k'); hold on;
semilogx([1000000 4500000],[-0.2 -0.2],'k'); hold on;
semilogx([6000000 10000000],[-51 -51],'k'); hold on;

semilogx(imag(p4)/(2*pi),a4,'k'); hold on;
text(4720000,-0.75,'IChebyA11'); hold on;
semilogx(imag(p5)/(2*pi),a5,'k'); hold on;
text(4420000,-1.5,'IChebyB12'); hold on;
hold off;

% Zobraznení nepropustného pasma 
figure(2)
semilogx(imag(p)/(2*pi),a1,'k');
text(4800000,-26,'CauerA6'); hold on;
axis([4500000 10000000 -80 0]);
xlabel('--> f[Hz]');
ylabel('--> a[dB]');
title('Zobrazení nepropustného pásma filtru');
hold on;
semilogx([4500000 4500000],[-80 -0.2],'k'); hold on;
semilogx([6000000 6000000],[-51 0],'k'); hold on;
semilogx([1000000 4500000],[-0.2 -0.2],'k'); hold on;
semilogx([6000000 10000000],[-51 -51],'k'); hold on;
semilogx(imag(p4)/(2*pi),a4,'k'); hold on;
text(7300000,-55,'IChebyA11'); hold on;
semilogx(imag(p5)/(2*pi),a5,'k'); hold on;
text(8100000,-77,'IChebyB12'); hold on;
hold off;

% Zobrazeni propustneho pasma 
figure(3)
semilogx(imag(p)/(2*pi),a,'k');
hold on;
text(4000000,-0.3,'CauerB6'); hold on;
axis([3000000 5000000 -2 1]);
xlabel('--> f[Hz]');
ylabel('--> a[dB]');
title('Zobrazení propustného pásma filtru');
hold on;
% Zobrazeni tolerancniho schematu
semilogx([4500000 4500000],[-80 -0.2],'k'); hold on;
semilogx([6000000 6000000],[-51 0],'k'); hold on;
semilogx([1000000 4500000],[-0.2 -0.2],'k'); hold on;
semilogx([6000000 10000000],[-51 -51],'k'); hold on;
semilogx(imag(p4)/(2*pi),a2,'k'); hold on;
text(4500000,0.1,'Cheby11'); hold on;
semilogx(imag(p5)/(2*pi),a3,'k'); hold on;
text(4700000,-1,'Butter26'); hold on;
hold off;

% Zobraznení nepropustného pasma 
figure(4)
semilogx(imag(p)/(2*pi),a,'k');
text(7500000,-55,'CauerB6'); hold on;
axis([4500000 10000000 -80 0]);
xlabel('--> f[Hz]');
ylabel('--> a[dB]');
title('Zobrazení nepropustného pásma filtru');
hold on;
semilogx([4500000 4500000],[-80 -0.2],'k'); hold on;
semilogx([6000000 6000000],[-51 0],'k'); hold on;
semilogx([1000000 4500000],[-0.2 -0.2],'k'); hold on;
semilogx([6000000 10000000],[-51 -51],'k'); hold on;
semilogx(imag(p4)/(2*pi),a2,'k'); hold on;
text(4550000,-26,'Cheby11'); hold on;
semilogx(imag(p5)/(2*pi),a3,'k'); hold on;
text(4850000,-3,'Butter26'); hold on;
hold off;



%---------------------------- Koeficienty filtru pro zobrazeni skupinového zpoždìní ------------------------

%----------------------------------------- Cauerova aproximace B6 ------------------------------------------
A = [1 41169680 2.159544E+15 5.410112E+22 1.210116E+30 1.579738E+37 1.258806E+44];
B = [29894420000000 0 1.255416E+29 0 1.230152E+44];
%----------------------------------------- Cauerova aproximace A6 ------------------------------------------
A1 = [1 41114920 2.174136E+15 5.471232E+22 1.238937E+30 1.6384E+37 1.358765E+44];
B1 = [0.002609761 0 48187470000000 0 1.552961E+29 0 1.327836E+44];
%--------------------------------------- Chebysevova aproximace 11 -----------------------------------------
A2 = [1 40674280 3.034456E+15 8.988043E+22 3.259484E+30 6.948896E+37 1.50746E+45 ...
      2.191208E+52 2.82191E+59 2.414811E+66 1.475948E+73 4.245405E+79];
B2 = [4.245405E+79];
%------------------- Skupinové zpozdìní pro Buttervothovu aproximaci nelze vypocitat -----------------------
 %A3 = [1 515605800 1.329247E+17 2.279386E+25 2.918159E+33 2.968156E+41 2.49244E+49 ...
 %     1.772788E+57 1.087373E+65 5.826349E+72 2.752842E+80 1.154677E+88 4.319795E+95 ...
 %     1.445669E+103 4.334089E+110 1.163992E+118 2.796438E+125 5.992165E+132 1.139874E+140 ...
 %     1.912039E+147 2.801772E+154 3.540451E+161 3.789469E+168 3.348563E+175 2.350851E+182 ...
 %     1.232171E+189 4.295768E+195 7.488259E+201];
 %B3 = [7.488257E+201];
%------------------------------------------ Invezní Chebysev A11 -------------------------------------------
A4 = [1 253849200 3.22195E+16 2.695297E+24 1.651391E+32 7.796332E+39 2.904446E+47 ...
      8.599991E+54 2.001274E+62 3.56485E+69 4.454311E+76 3.420419E+83];
B4 = [646478.8 0 1.830239E+22 0 1.45084E+38 0 4.694236E+53 0 6.6449E+68 0 3.420419E+83];
%------------------------------------------ Inverzní Chebysev B12 ------------------------------------------
A5 = [1 273404500 3.7375E+16 3.373361E+24 2.238477E+32 1.151658E+40 4.722024E+47 ...
      1.561028E+55 4.157556E+62 8.782616E+69 1.423261E+77 1.623672E+84 1.102774E+91];
B5 = [46110810000000 0 9.320985E+29 0 6.118703E+45 0 1.752887E+61 0 2.279522E+76 0 ...
      1.102774E+91];
%------------------------------------------- Invezní Chebysev A13 ------------------------------------------
A11 = [1 302751600 4.582893E+16 4.588844E+24 3.391115E+32 1.95538E+40 9.069643E+47 ...
       3.439874E+55 1.071819E+63 2.736478E+70 5.611579E+77 9.016487E+84 1.027293E+92 ...
       7.59185E+98];
B11 = [830655.5 0 2.777563E+22 0 2.653621E+38 0 1.086518E+54 0 2.162572E+69 0 ...
       2.066072E+84 0 7.591839E+98];
%------------------------------------------ Cauerova aproximace B8 ----------------------------------------- 
A12 = [1 40615010 2.624872E+15 7.321153E+22 2.205082E+30 4.031134E+37 6.447161E+44 ...
       6.356483E+51 3.771198E+58];
B12 = [6153305000000 0 3.622615E+28 0 6.515352E+43 0 3.685356E+58];
%---------------------------------------- Cauerova aproximace B8 5MHz --------------------------------------
A13 = [1 40495040 2.715023E+15 7.684398E+22 2.411392E+30 4.560935E+37 7.727922E+44 ...
       8.132186E+51 5.537859E+58];
B13 = [19916260000000 0 8.824821E+28 0 1.223489E+44 0 5.411802E+58];
%---------------------------- Nejde zobrazit skupinové spozdeni pro ICheby A18 -----------------------------
 %A14 = [1 431251200 9.298946E+16 1.331628E+25 1.419211E+33 1.195809E+41 8.26235E+48 ...
 %       4.791712E+56 2.369116E+64 1.007989E+72 3.712994E+79 1.184697E+87 3.273525E+94 ...
 %       7.754676E+101 1.566086E+109 2.605239E+116 3.538005E+123 3.42178E+130 2.392108E+137];
 %B14 = [0.002053873 0 3.283888E+14 0 8.642845E+30 0 8.757615E+46 0 4.445187E+62 ...
 %       0 1.26742E+78 0 2.122728E+93 0 2.072026E+108 0 1.09067E+123 0 2.392108E+137]; 

delay = grdelays(B,A,imag(p));     		   % Cauer B6	6MHz   		
delay2 = grdelays(B2,A2,imag(p2));			% Chebysev 11
delay4 = grdelays(B4,A4,imag(p4));			% ICheby A11
delay5 = grdelays(B5,A5,imag(p5));			% ICheby B12
	

%delay11 = grdelays(B11,A11,imag(p11));	% ICheby A13 5.5MHz
%delay12 = grdelays(B12,A12,imag(p12));	% Cauer B8 5.5MHz
%delay13 = grdelays(B13,A13,imag(p13));	% Cauer B8 5MHz	
%delay14 = grdelays(B14,A14,imag(p14));	% ICheby A18 5MHz


% Zobrazení celeho skupinoveho zpozdeni 

figure(5);
plot(imag(p)/(2*pi),delay,'k');
xlabel('--> f[Hz]');
ylabel('--> t[s]');
hold on;
plot([4500000 4500000],[0 1e-6],'k');
hold on;
title('Zobrazení propustného pásma skupinového zpoždìní');
hold on;
plot(imag(p)/(2*pi),delay,'k');
text(5.4e+6,0.15e-6,'CauerB6'); hold on;
axis([0 10000000 0 10e-7]);

plot(imag(p2)/(2*pi),delay2,'k');hold on;
text(3.2e+6,0.95e-6,'Cheby11'); hold on;
plot(imag(p4)/(2*pi),delay4,'k');hold on;
plot(imag(p5)/(2*pi),delay5,'k');hold on;
text(5.2e+6,0.65e-6,'IChebyB12/A11'); hold off;
grid; 



figure(6);
plot(imag(p)/(2*pi),delay,'k');
xlabel('--> f[Hz]');
ylabel('--> t[s]');
hold on;
plot([4500000 4500000],[0 1e-6],'k');
hold on;
title('Zobrazení propustného pásma skupinového zpoždìní');
hold on;

plot(imag(p)/(2*pi),delay,'k');
text(3.9e+6,5.25e-7,'CauerB6'); hold on;
axis([0 4500000 1e-7 6e-7]);
plot(imag(p4)/(2*pi),delay4,'k');hold on;
text(3.7e+6,2.25e-7,'IChebyA11'); 
plot(imag(p5)/(2*pi),delay5,'k');hold on;
text(0.25e+6,1.75e-7,'IChebyB12'); 
grid;



