% Example sensitivity
% Course: Analog and Digital Filter Design
% Lecture: Specification of an analog filter
% Authors: Miroslav D. Lutovac, Dejan V. Tosic
% Email: lutovac@kondor.etf.bg.ac.yu
% Web: http://kondor.etf.bg.ac.yu/~lutovac/
% Copyright (c) 2005 by Lutovac & Tosic
% $Revision: 1.00 $  $Date: 2005/01/30$
% See also: Miroslav D. Lutovac, Dejan V. Tosic, Brian L. Evans
%  Filter Design for Signal Processing Using MATLAB and Mathematica
%  Prentice Hall (c)2001. ISBN 0-201-36130-2

clear all, close all, clc

Fp = 3400;
Fs = 5000;
Ap = 1;
As = 20;

[n, Wn] = cheb1ord(Fp, Fs, Ap, As, 's')
[z,p,k] = cheb1ap(n, Ap);
Q = (-abs(p)./(2*real(p)));
w = round(abs(p)*(2*pi*Fp));
D1 = [1 w(1)/Q(1) w(1)^2];
D2 = [1 w(2)/Q(2) w(2)^2];
N1 = w(1)^2;
N2 = w(2)^2;
D = conv(D1,D2);
N = 10^(-Ap/20)*conv(N1,N2);

Fplot = 6000;
f = 0:Fplot/200:Fplot;
H0 = freqs(N,D,2*pi*f);
M0 = -20*log10(abs(H0));

tolerance = 10;
Q = (1 + tolerance/100)*(-abs(p)./(2*real(p)));
w = (1 + tolerance/100)*round(abs(p)*(2*pi*Fp));
D1 = [1 w(1)/Q(1) w(1)^2];
D2 = [1 w(2)/Q(2) w(2)^2];
N1 = w(1)^2;
N2 = w(2)^2;
D = conv(D1,D2);
N = 10^(-Ap/20)*conv(N1,N2);
H1 = freqs(N,D,2*pi*f);
M1 = -20*log10(abs(H1));

Q = (1 - tolerance/100)*(-abs(p)./(2*real(p)));
w = (1 - tolerance/100)*round(abs(p)*(2*pi*Fp));
D1 = [1 w(1)/Q(1) w(1)^2];
D2 = [1 w(2)/Q(2) w(2)^2];
N1 = w(1)^2;
N2 = w(2)^2;
D = conv(D1,D2);
N = 10^(-Ap/20)*conv(N1,N2);
H2 = freqs(N,D,2*pi*f);
M2 = -20*log10(abs(H2));

plot(f,M0,f,M1,'--',f,M2,':',[0 Fp Fp],[Ap Ap As],'r--',[Fs Fs Fplot],[0 As As],'r--')
xlabel('f (Hz)')
ylabel('A (dB)')
legend('ideal','+10%','-10%')
