# ellip
椭圆滤波器的设计
函数库: TySignalProcessing
注意
不同硬件环境下，由于计算精度的问题导致运行结果可能不一致。
# 语法
b, a = ellip(n, Rp, Rs, Wp)
b, a = ellip(n, Rp, Rs, Wp, ftype)
z, p, k = ellip(_____; otype = "zpk")
A, B, C, D = ellip(___; otype = "ABCD")
___ = ellip(___, "s"; otype = otype)
# 说明
b, a = ellip(n, Rp, Rs, Wp) 返回具有归一化通带边缘频率 Wp 和峰峰值通带纹波的 Rp 分贝的 n 阶低通椭圆滤波器的传递函数系数。
b, a = ellip(n, Rp, Rs, Wp, ftype) 根据 ftype 的值和 Wp 的元素数量，设计低通、高通、带通或带阻椭圆滤波器。由此产生的带通和带阻设计为 2n 阶。
z, p, k = ellip(_____; otype = "zpk") 设计一个低通、高通、带通或带阻数字椭圆滤波器并返回其零点、极点和增益。这个语法可以包括前面语法中的任何输入参数。
A, B, C, D = ellip(___; otype = "ABCD") 设计一个低通、高通、带通或带阻数字椭圆滤波器并返回指定其状态空间表示的矩阵。
___ = ellip(___, "s"; otype = otype) 设计一个低通、高通、带通或带阻的模拟椭圆滤波器，其通带边缘角频率为 Wp，通带纹波为 Rp 分贝，阻带衰减为 Rs 分贝。
# 示例
低通椭圆传递函数
设计一个具有 5 dB 通带纹波、40 dB 阻带衰减和 300 Hz 通带边缘频率的 6 阶低通椭圆滤波器，对于以 1000 Hz 采样的数据，对应于 0.6π rad/sample。绘制其幅度和相位响应。用它来对 1000 个样本的随机信号进行滤波。
using TyPlot
using TyMath
using TySignalProcessing
fc = 300
fs = 1000
b, a = ellip(6, 10, 50, fc / (fs / 2))
freqz(b, a, [], fs; plotfig=true)
subplot(2, 1, 1)
ylim([-100 20])
rng = MT19937ar(1234)
dataIn = randn(rng, 1000, 1)
dataOut, = filter1(b, a, dataIn)
带阻椭圆滤波器
设计一个 6 阶椭圆带阻滤波器，其归一化边缘频率为 0.2π 和 0.6π rad/sample、5 dB 通带纹波和 50 dB 阻带衰减。绘制其幅度和相位响应。用它来对随机数据进行滤波。
using TyMath
using TySignalProcessing
b, a = ellip(3, 5, 50, [0.2 0.6], "bandstop")
freqz(b, a; plotfig=true)
rng = MT19937ar(1234)
dataIn = randn(rng, 1000, 1)
dataOut, = filter1(b, a, dataIn)
高通椭圆滤波器
设计一个通带边缘频率为 300 Hz 的 6 阶高通椭圆滤波器，对于以 1000 Hz 采样的数据，对应于 0.6π rad/sample。指定 3 dB 的通带纹波和 50 dB 的阻带衰减。绘制幅度和相位响应。
using TySignalProcessing
z, p, k = ellip(6, 3, 50, 300 / 500, "highpass"; otype="zpk")
sos = zp2sos(z, p, k; g_flag=false)
fvtool(sos; Analysis="freq")
带通椭圆滤波器
设计一个椭圆带通滤波器，其较低的通带频率为 500 Hz，较高的通带频率为 560 Hz。指定 3 dB 的通带纹波、40 dB 的阻带衰减。
using TySignalProcessing
b, a = ellip(10, 3, 40, [500 560] / 750, "bandpass")
fvtool(b, a)
模拟 IIR 低通滤波器对比
设计一个截止频率为 2 GHz 的 5 阶模拟巴特沃斯低通滤波器。乘以 2π，将频率转换为 rad/s。计算 4096 点处滤波器的频率响应。
using TyPlot
using TySignalProcessing
using TyControlSystems
n = 5
f = 2e9
zb, pb, kb = butter(n, 2 * pi * f, "s"; otype="zpk")
bb, ab = zp2tf(zb, pb, kb)
hb, wb = freqs(bb, ab, 4096)
设计具有相同边缘频率和 3dB 通带纹波的5阶切比雪夫 I 型滤波器。计算其频率响应。
z1, p1, k1 = cheby1(n, 3, 2 * pi * f, "s"; otype="zpk")
b1, a1 = zp2tf(z1, p1, k1)
h1, w1 = freqs(b1, a1, 4096)
设计具有相同边缘频率和 30dB 阻带衰减的5阶切比雪夫 II 型滤波器。计算其频率响应。
z2, p2, k2 = cheby2(n, 30, 2 * pi * f, "s"; otype="zpk")
b2, a2 = zp2tf(z2, p2, k2)
h2, w2 = freqs(b2, a2, 4096)
设计一个边缘频率相同、通带纹波为 3dB、阻带衰减为 30dB 的 5 阶椭圆滤波器。计算其频率响应。
ze, pe, ke = ellip(n, 3, 30, 2 * pi * f, "s"; otype="zpk")
be, ae = zp2tf(ze, pe, ke)
he, we = freqs(be, ae, 4096)
设计一个具有相同边缘频率的 5 阶贝塞尔滤波器。计算其频率响应。
zf, pf, kf = besself(n, 2 * pi * f, "lowpass", otype = "zpk")
bf, af = zp2tf(zf, pf, kf)
hf, wf = freqs(bf, af, 4096)
以分贝为单位绘制衰减曲线。以千兆赫表示频率。比较滤波器。
plot(wb / (2e9 * pi), mag2db.(abs.(hb)))
hold("on")
plot(w1 / (2e9 * pi), mag2db.(abs.(h1)))
plot(w2 / (2e9 * pi), mag2db.(abs.(h2)))
plot(we / (2e9 * pi), mag2db.(abs.(he)))
plot(wf / (2e9 * pi), mag2db.(abs.(hf)))
axis([0 5 -45 5])
grid("on")
xlabel("Frequency (GHz)")
ylabel("Attenuation (dB)")
legend(["butter", "cheby1", "cheby2", "ellip", "besself"])
巴特沃斯和切比雪夫 II 型滤波器的通带平坦，过渡带宽。切比雪夫 I 型和椭圆型滤波器的滚降速度较快，但有通带纹波。切比雪夫 II 型设计函数的频率输入是设定阻带的起始点，而不是通带的终点。贝塞尔滤波器在通带上有近似恒定的群延迟。
# 输入参数
n - 滤波器阶数整数标量
滤波器阶数，指定为整数标量。对于带通和带阻设计，n 表示滤波器阶数的一半。
数据类型： Int
Rp - 峰值到峰值通带纹波正标量
峰值到峰值通带纹波，指定为在分贝中表示的正标量。
数据类型： Float
Rs - 阻带衰减正标量
从通带峰值向下的阻带衰减，指定为以分贝表示的正标量。
如果您的规格，ℓ, 是线性单位，可以使用
数据类型： Float
Wp - 通带边缘频率标量 | 两个元素的向量
通带边缘频率，指定为标量或双元素向量。通带边沿频率是滤波器幅度响应为 -RP 分贝的频率。较小的通带纹波，RP 值，导致更宽的过渡带。 当为两个元素的向量时，需要指定滤波器类型：bandpass，bandstop
数据类型： Float
ftype - 滤波器类型"lowpass"(默认) | "bandpass" | "highpass" | "bandstop"
滤波器的类型。默认为 "lowpass"
数据类型： String
otype - 输出类型"ba"(默认) | "zpk" | "ABCD"
滤波器表示类型。默认为 "ba"。
数据类型： String
# 输出参数
b, a - 传递函数系数向量
IIR 滤波器的分子（B）和分母（A）多项式。
数据类型： Float
z, p, k - 零点，极点和增益向量，标量
零点，极点和滤波器的增益，返回为长度 N 的两个向量（带通和带阻滤波器时，长度为 2N）和一个标量。
数据类型： Float
A, B, C, D - 状态空间矩阵矩阵
滤波器的状态空间表示，以矩阵形式返回。如果低通和高通设计 m = n，带通和带阻滤波器 m = 2n，则 A 为 m × m, B 为 m × 1, C 为 1 × m, D 为 1 × 1。
对于数字滤波器，状态空间矩阵将状态向量 x、输入 u 和输出 y 通过以下式子联系起来：
对于模拟滤波器，状态空间矩阵将状态向量 x、输入 u 和输出 y 通过以下式子联系起来：
数据类型： Float
# 详细信息
椭圆滤波器提供比巴特沃斯或切比雪夫滤波器更陡峭的滚降特性，但在通带和阻带中都是等波纹的。通常，椭圆滤波器满足给定的性能规格，并且是任何滤波器类型中最低阶的。
ellip 使用五步算法：
它使用 ellipap 来找到低通模拟滤波器的模板零点，极点和增益；
它将零点，极点和增益转换为状态空间形式；
如果需要，它使用状态空间变换，将低通滤波器转换为具有所需频率约束的高通，带通或带阻滤波器；
对于数字滤波器设计，它使用双线性将模拟滤波器转换为数字滤波器，通过双线性转换与频率预调高。仔细的频率调节使模拟滤波器和数字滤波器能够在 WP 或 W1 和 W2 处具有相同的频率响应幅度；
它根据需要将状态空间滤波器转换回传输函数或零极点形式。
# 版本历史记录
在 2024a 之前推出
2025a SP3：零极点输出统一为复数向量
z,p,k=ellip(__;otype="zpk") 用法，输出零点 z 由单列矩阵或向量修改为始终为复数向量。
fc = 300
fs = 1000
z,p,k = ellip(6, 10, 50, fc / (fs / 2);otype="zpk")
更新前结果：
z
6×1 Matrix{ComplexF64}:
-0.8614632003131223 + 0.507820002073837im
-0.8614632003131223 - 0.507820002073837im
-0.3710290323059894 + 0.9286212668176853im
-0.3710290323059894 - 0.9286212668176853im
-0.48055349687223803 + 0.876965413596148im
-0.48055349687223803 - 0.876965413596148im
2025a SP3 更新后结果：
z
6-element Vector{ComplexF64}:
-0.861463200313123 + 0.5078200020738369im
-0.861463200313123 - 0.5078200020738369im
-0.3710290323059896 + 0.9286212668176846im
-0.3710290323059896 - 0.9286212668176846im
-0.4805534968722369 + 0.8769654135961474im
-0.4805534968722369 - 0.8769654135961474im