1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55
|
// Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
// Copyright (C) INRIA - 1989 - F.Delebecque
// Copyright (C) INRIA - 1996 - C. Bunks
//
// This file must be used under the terms of the CeCILL.
// This source file is licensed as described in the file COPYING, which
// you should have received as part of this distribution. The terms
// are also available at
// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
function [ze,po,gain]=zpell(epsilon,A,omegac,omegar)
//[ze,po,gain]=zpell(epsilon,A,omegac,omegar)
//Poles and zeros of prototype lowpass elliptic filter
//gain is the gain of the filter
// epsilon :Ripple of filter in pass band (0<epsilon<1)
// A :Attenuation of filter in stop band (A>1)
// omegac :Pass band cut-off frequency in Hertz
// omegar :Stop band cut-off frequency in Hertz
// ze :Resulting zeros of filter
// po :Resulting poles of filter
// gain :Resulting gain of filter
//
//!
m1=(epsilon*epsilon)/(A*A-1);
K1=%asn(1,m1);
K1t=imag(%asn(1/sqrt(m1),m1));
m=(omegac/omegar)^2;
K=%asn(1,m);
Kt=imag(%asn(1/sqrt(m),m));
n=(K1t*K)/(K1*Kt);
order=round(n);
u0=-(Kt/K1t)*%asn(sqrt(1/(1+epsilon*epsilon)),1-m1);
even=2*int(order/2);
if order<>even then,
vmin=2*K/n;
else,
vmin=K/n;
end,
v=vmin:(2*K/n):K;
un=ones(1:max(size(v)));
zlambda=-un*Kt+%i*v;
plambda= u0*un+%i*v;
ze=%i*imag(%i*omegac*%sn(-%i*zlambda,m));
ze=[ze,conj(ze)];
po=%i*omegac*%sn(-%i*plambda,m);
po=[po,conj(po)];
if order<>even then,
po=[po,%i*omegac*%sn(-%i*u0,m)];
end,
gain=abs(real(prod(po))/real(prod(ze)));
if order==even then,
gain=gain/sqrt(1+epsilon^2);
end;
endfunction
|