File: zpell.sci

package info (click to toggle)
scilab 5.3.3-10
  • links: PTS, VCS
  • area: main
  • in suites: wheezy
  • size: 330,656 kB
file content (55 lines) | stat: -rw-r--r-- 1,731 bytes parent folder | download
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