零阶保持器（ZOH）变换和Tustin离散化变换以及可视化
引言
-
零阶保持器（ZOH）变换：
基于在采样间隔内保持采样值不变的思想。在每个采样周期 T s T_s Ts 内，将连续信号的值保持为该周期开始时刻的采样值。从数学角度，先将连续传递函数 G ( s ) G(s) G(s) 与零阶保持器 G h ( s ) = 1 − e − T s s s G_h(s) = \frac{1 - e^{-T_s s}}{s} Gh(s)=s1−e−Tss 串联，再对 G h ( s ) G ( s ) G_h(s)G(s) Gh(s)G(s) 进行 z z z 变换，即
G ( z ) = ( 1 − z − 1 ) Z [ G ( s ) s ] G(z) = (1 - z^{-1})\mathcal{Z}\left[\frac{G(s)}{s}\right] G(z)=(1−z−1)Z[sG(s)] -
Tustin离散化：
基于对拉普拉斯域 s s s 与 z z z 域之间的一种映射关系，用梯形积分法近似连续时间系统的积分环节，其变换公式为
s = 2 T s 1 − z − 1 1 + z − 1 s = \frac{2}{T_s}\frac{1 - z^{-1}}{1 + z^{-1}} s=Ts21+z−11−z−1
将连续系统传递函数 G ( s ) G(s) G(s) 中的 s s s 用此式替换来得到离散系统传递函数 G ( z ) G(z) G(z)。
下面将详细介绍零阶保持器（ZOH）变换的步骤，并通过一个具体例子来说明如何进行变换。
零阶保持器（ZOH）变换
零阶保持器（ZOH）变换用于将连续时间系统转换为离散时间系统，其主要步骤如下：
- 确定零阶保持器的传递函数：零阶保持器的传递函数为 G Z O H ( s ) = 1 − e − s T s G_{ZOH}(s)=\frac{1 - e^{-sT}}{s} GZOH(s)=s1−e−sT，其中 T T T 是采样周期。
- 将零阶保持器与连续系统传递函数相乘：设连续系统的传递函数为 G ( s ) G(s) G(s)，那么相乘后的结果为 G e q ( s ) = G Z O H ( s ) G ( s ) = 1 − e − s T s G ( s ) G_{eq}(s)=G_{ZOH}(s)G(s)=\frac{1 - e^{-sT}}{s}G(s) Geq(s)=GZOH(s)G(s)=s1−e−sTG(s)。
- 对 G e q ( s ) G_{eq}(s) Geq(s) 进行 Z Z Z 变换：得到离散时间系统的传递函数 G ( z ) = Z [ G e q ( s ) ] G(z)=Z[G_{eq}(s)] G(z)=Z[Geq(s)]。
具体例子
假设我们有一个简单的连续时间系统，其传递函数为 G ( s ) = 1 s + a G(s)=\frac{1}{s + a} G(s)=s+a1（ a a a 为常数），采样周期为 T T T，下面我们来对其进行 ZOH 变换。
步骤 1：确定零阶保持器的传递函数
零阶保持器的传递函数 G Z O H ( s ) = 1 − e − s T s G_{ZOH}(s)=\frac{1 - e^{-sT}}{s} GZOH(s)=s1−e−sT。
步骤 2：将零阶保持器与连续系统传递函数相乘
G e q ( s ) = G Z O H ( s ) G ( s ) = 1 − e − s T s ⋅ 1 s + a = ( 1 − e − s T ) 1 s ( s + a ) G_{eq}(s)=G_{ZOH}(s)G(s)=\frac{1 - e^{-sT}}{s}\cdot\frac{1}{s + a}=(1 - e^{-sT})\frac{1}{s(s + a)} Geq(s)=GZOH(s)G(s)=s1−e−sT⋅s+a1=(1−e−sT)s(s+a)1
我们先对
1
s
(
s
+
a
)
\frac{1}{s(s + a)}
s(s+a)1 进行部分分式展开：
设
1
s
(
s
+
a
)
=
A
s
+
B
s
+
a
\frac{1}{s(s + a)}=\frac{A}{s}+\frac{B}{s + a}
s(s+a)1=sA+s+aB
通分可得：
1
=
A
(
s
+
a
)
+
B
s
=
(
A
+
B
)
s
+
A
a
1 = A(s + a)+Bs=(A + B)s+Aa
1=A(s+a)+Bs=(A+B)s+Aa
通过系数相等可得：
{
A
+
B
=
0
A
a
=
1
\begin{cases}A + B = 0\\Aa = 1\end{cases}
{A+B=0Aa=1
解得
A
=
1
a
A=\frac{1}{a}
A=a1，
B
=
−
1
a
B=-\frac{1}{a}
B=−a1
所以
1
s
(
s
+
a
)
=
1
a
(
1
s
−
1
s
+
a
)
\frac{1}{s(s + a)}=\frac{1}{a}(\frac{1}{s}-\frac{1}{s + a})
s(s+a)1=a1(s1−s+a1)
则 G e q ( s ) = ( 1 − e − s T ) 1 a ( 1 s − 1 s + a ) = 1 a ( 1 s − 1 s + a ) − 1 a e − s T ( 1 s − 1 s + a ) G_{eq}(s)=(1 - e^{-sT})\frac{1}{a}(\frac{1}{s}-\frac{1}{s + a})=\frac{1}{a}(\frac{1}{s}-\frac{1}{s + a})-\frac{1}{a}e^{-sT}(\frac{1}{s}-\frac{1}{s + a}) Geq(s)=(1−e−sT)a1(s1−s+a1)=a1(s1−s+a1)−a1e−sT(s1−s+a1)
步骤 3：对 G e q ( s ) G_{eq}(s) Geq(s) 进行 Z Z Z 变换
根据 Z Z Z 变换的性质和常用 Z Z Z 变换对：
- Z [ 1 s ] = z z − 1 Z[\frac{1}{s}]=\frac{z}{z - 1} Z[s1]=z−1z
- Z [ 1 s + a ] = z z − e − a T Z[\frac{1}{s + a}]=\frac{z}{z - e^{-aT}} Z[s+a1]=z−e−aTz
- 以及时移定理 Z [ e − s T F ( s ) ] = z − 1 Z [ F ( s ) ] Z[e^{-sT}F(s)]=z^{-1}Z[F(s)] Z[e−sTF(s)]=z−1Z[F(s)]
对于
1
a
(
1
s
−
1
s
+
a
)
\frac{1}{a}(\frac{1}{s}-\frac{1}{s + a})
a1(s1−s+a1)，其
Z
Z
Z 变换为：
1
a
(
z
z
−
1
−
z
z
−
e
−
a
T
)
=
z
a
(
z
−
e
−
a
T
)
−
(
z
−
1
)
(
z
−
1
)
(
z
−
e
−
a
T
)
=
z
(
1
−
e
−
a
T
)
a
(
z
−
1
)
(
z
−
e
−
a
T
)
\frac{1}{a}(\frac{z}{z - 1}-\frac{z}{z - e^{-aT}})=\frac{z}{a}\frac{(z - e^{-aT})-(z - 1)}{(z - 1)(z - e^{-aT})}=\frac{z(1 - e^{-aT})}{a(z - 1)(z - e^{-aT})}
a1(z−1z−z−e−aTz)=az(z−1)(z−e−aT)(z−e−aT)−(z−1)=a(z−1)(z−e−aT)z(1−e−aT)
对于
1
a
e
−
s
T
(
1
s
−
1
s
+
a
)
\frac{1}{a}e^{-sT}(\frac{1}{s}-\frac{1}{s + a})
a1e−sT(s1−s+a1)，其
Z
Z
Z 变换为：
z
−
1
z
(
1
−
e
−
a
T
)
a
(
z
−
1
)
(
z
−
e
−
a
T
)
=
1
−
e
−
a
T
a
(
z
−
1
)
(
z
−
e
−
a
T
)
z^{-1}\frac{z(1 - e^{-aT})}{a(z - 1)(z - e^{-aT})}=\frac{1 - e^{-aT}}{a(z - 1)(z - e^{-aT})}
z−1a(z−1)(z−e−aT)z(1−e−aT)=a(z−1)(z−e−aT)1−e−aT
那么 G ( z ) = z ( 1 − e − a T ) a ( z − 1 ) ( z − e − a T ) − 1 − e − a T a ( z − 1 ) ( z − e − a T ) = ( z − 1 ) ( 1 − e − a T ) a ( z − 1 ) ( z − e − a T ) = 1 − e − a T a ( z − e − a T ) G(z)=\frac{z(1 - e^{-aT})}{a(z - 1)(z - e^{-aT})}-\frac{1 - e^{-aT}}{a(z - 1)(z - e^{-aT})}=\frac{(z - 1)(1 - e^{-aT})}{a(z - 1)(z - e^{-aT})}=\frac{1 - e^{-aT}}{a(z - e^{-aT})} G(z)=a(z−1)(z−e−aT)z(1−e−aT)−a(z−1)(z−e−aT)1−e−aT=a(z−1)(z−e−aT)(z−1)(1−e−aT)=a(z−e−aT)1−e−aT
这样我们就完成了对连续系统 G ( s ) = 1 s + a G(s)=\frac{1}{s + a} G(s)=s+a1 的 ZOH 变换，得到了离散时间系统的传递函数 G ( z ) = 1 − e − a T a ( z − e − a T ) G(z)=\frac{1 - e^{-aT}}{a(z - e^{-aT})} G(z)=a(z−e−aT)1−e−aT。
Tustion离散化
Tustin离散化，又称双线性变换，是一种将连续时间系统转换为离散时间系统的方法，在控制系统设计等领域应用广泛，以下从原理、公式及特点进行介绍：
原理
Tustin离散化基于对拉普拉斯域 s s s 与 z z z 域之间的一种映射关系，通过对连续时间系统的传递函数进行变换，得到对应的离散时间系统传递函数。从本质上来说，它是用梯形积分法来近似连续时间系统的积分环节，从而实现从连续到离散的转换。
变换公式
s
s
s 域与
z
z
z 域的映射关系为：
s
=
2
T
s
1
−
z
−
1
1
+
z
−
1
s = \frac{2}{T_s}\frac{1 - z^{-1}}{1 + z^{-1}}
s=Ts21+z−11−z−1
其中 T s T_s Ts 是采样周期。如果已知连续时间系统的传递函数 G ( s ) G(s) G(s)，将其中的 s s s 用上述表达式替换，经过化简后就可以得到离散时间系统的传递函数 G ( z ) G(z) G(z)。
Tustin变换特点
-
优点：
- Tustin变换是一种保角映射，它能将 s s s 平面的左半平面映射到 z z z 平面的单位圆内，因此如果连续系统是稳定的，经过Tustin变换得到的离散系统也通常是稳定的，这对于控制系统的稳定性保证很重要。
- 它对于具有有理函数形式的连续传递函数，能得到精确的有理形式的离散传递函数。
-
缺点：
- 存在频率畸变现象，它会对高频段的频率响应产生较大的失真。为了补偿这种畸变，在使用Tustin变换时，常常需要进行预畸变处理，对截止频率等参数进行调整。
以下是常见的双线性变换（Tustin变换）结论：
1. 变换关系
连续域
s
s
s 与离散域
z
z
z 之间的映射关系为
s
=
2
T
s
1
−
z
−
1
1
+
z
−
1
s = \frac{2}{T_s}\frac{1 - z^{-1}}{1 + z^{-1}}
s=Ts21+z−11−z−1
其中
T
s
T_s
Ts 是采样周期。
反之，也可由
z
=
1
+
T
s
2
s
1
−
T
s
2
s
z = \frac{1 + \frac{T_s}{2}s}{1 - \frac{T_s}{2}s}
z=1−2Tss1+2Tss
从
s
s
s 得到
z
z
z。
2. 频率变换特性
-
低频特性：
在低频段，Tustin变换的频率映射近似线性，能够较好地保持连续系统的频率响应特性。 -
高频特性：
存在频率畸变现象，会将连续系统的无限大频率范围映射到离散系统的 ± π T s \pm\frac{\pi}{T_s} ±Tsπ （数字频率 ± π \pm\pi ±π）范围内，导致高频段出现失真。为补偿畸变，常进行预畸变处理，即对连续系统截止频率等参数预先调整。
3. 稳定性
Tustin变换将
s
s
s 平面的左半平面（
Re
(
s
)
<
0
\text{Re}(s)<0
Re(s)<0）映射到
z
z
z 平面的单位圆内（
∣
z
∣
<
1
\vert z\vert<1
∣z∣<1），
右半平面（
Re
(
s
)
>
0
\text{Re}(s)>0
Re(s)>0）映射到单位圆外（
∣
z
∣
>
1
\vert z\vert>1
∣z∣>1），
虚轴（
Re
(
s
)
=
0
\text{Re}(s)=0
Re(s)=0）映射到单位圆上（
∣
z
∣
=
1
\vert z\vert = 1
∣z∣=1）。
因此，若连续系统稳定（极点在 s s s 左半平面），经Tustin变换后的离散系统也稳定（极点在 z z z 平面单位圆内）。
4. 积分器变换
连续积分器
G
(
s
)
=
1
s
G(s)=\frac{1}{s}
G(s)=s1 经Tustin变换后，离散传递函数为
G
(
z
)
=
T
s
2
z
+
1
z
−
1
G(z) = \frac{T_s}{2}\frac{z + 1}{z - 1}
G(z)=2Tsz−1z+1
这种变换是基于用梯形积分近似连续时间积分环节得到的。
5. 微分器变换
连续微分器
G
(
s
)
=
s
G(s)=s
G(s)=s 经Tustin变换后，离散传递函数为
G
(
z
)
=
2
T
s
z
−
1
z
+
1
G(z) = \frac{2}{T_s}\frac{z - 1}{z + 1}
G(z)=Ts2z+1z−1
在Matlab中，可使用c2d函数将连续时间域（s域）的系统转换为离散时间域（z域）的系统，从而实现s域传递函数到z域传递函数的变换。其语法为：
sysd = c2d(sysc, Ts, method)
其中：
sysc：是连续时间域的系统模型，一般通过tf函数创建。例如，若有连续传递函数 G ( s ) = 400 s 2 + 50 s G(s) = \frac{400}{s^2 + 50s} G(s)=s2+50s400，可使用sysc = tf(400,[1,50,0])来创建其模型。Ts：是采样时间，为一个正实数，表示采样周期。method：是离散化方法，常用取值如下：'zoh'：零阶保持（Zero - Order Hold）法，在每个采样间隔内，将连续信号的值保持为该周期开始时刻的采样值。'foh'：一阶保持（First - Order Hold）法，利用线性外推来近似采样间隔内的信号。'tustin'或'bilinear'：双线性变换（Tustin变换）法，基于拉普拉斯域 s s s与 z z z域之间的映射关系
s = 2 T s 1 − z − 1 1 + z − 1 s = \frac{2}{T_s}\frac{1 - z^{-1}}{1 + z^{-1}} s=Ts21+z−11−z−1
来进行离散化。
示例代码
% 创建连续传递函数模型
sysc = tf(1, [1, 1]); % G(s) = 1/(s + 1)
Ts = 0.1; % 采样周期为0.1
% 离散化：零阶保持法与双线性变换法
sysd_zoh = c2d(sysc, Ts, 'zoh');
sysd_tustin = c2d(sysc, Ts, 'tustin');
% 定义频率范围
freq_range = {1e-2, 1e2};
% 创建图窗并绘制伯德图
figure;
hold on;
bodeplot(sysc, 'b', freq_range); % 连续系统
bodeplot(sysd_zoh, 'r--', freq_range); % 零阶保持法离散化
bodeplot(sysd_tustin, 'g-.', freq_range); % 双线性变换法离散化
grid on;
legend('连续系统 G(s)', '零阶保持法离散化 G(z)', '双线性变换法离散化 G(z)', 'Location', 'Best');
title('连续与离散系统伯德图频率响应对比');
% 设置线宽（对所有线对象设置）
set(findall(gcf, 'Type', 'Line'), 'LineWidth', 1.5);
% 正确设置坐标轴字体大小：获取所有的坐标轴对象
ax = findobj(gcf, 'Type', 'axes');
set(ax, 'FontSize', 12);
% 设置图窗背景及大小
set(gcf, 'Color', 'w', 'Position', [100, 100, 800, 600]);
可以观察到以下区别：
幅值响应区别
- 低频段：连续系统与两种离散化方法（零阶保持法和双线性变换法）得到的离散系统幅值响应较为接近。这表明在低频段，两种离散化方法都能较好地保留连续系统的幅值特性，系统增益变化趋势相似。
- 高频段：零阶保持法离散化的系统幅值响应在高频段下降速度较快，呈现出低通滤波特性。相比之下，双线性变换法离散化的系统幅值响应在高频段下降相对较缓，并且在一定频率后出现了较为明显的差异。这体现了双线性变换法的频率畸变特性，它将连续系统的高频范围映射到离散系统有限的频率范围内，导致高频段幅值响应与连续系统有所不同。
相位响应区别
- 低频段：连续系统和离散系统的相位响应在低频段也较为相似，相位变化趋势相近，说明低频段相位特性离散化后变化不大。
- 高频段：零阶保持法离散化的系统相位滞后比连续系统更为严重，并且随着频率升高，相位下降速度加快。双线性变换法离散化的系统相位响应在高频段与连续系统差异也逐渐增大，但变化趋势与零阶保持法有所不同，其相位变化相对较为平缓一些，但同样存在明显的相位偏差。
总体而言，两种离散化方法在低频段对连续系统频率响应的近似效果较好，而在高频段，由于各自离散化原理的不同，与连续系统的频率响应均产生了不同程度的偏差。