自动控制一阶线性系统参数测定

本文链接：https://blog.csdn.net/qqttd281418216/article/details/120279550

本文介绍了如何测定一阶线性系统的参数，用于PID温控电烙铁的设计。通过微分方程和傅里叶变换分析系统，利用最小二乘法和线性近似处理噪声数据，推导出G和C的计算公式。同时，讨论了参数选择对精度的影响以及使用加权平均估计的注意事项。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

想做一个PID温控电烙铁，首先需要测定系统参数，便于以后的仿真和选择PID参数
设电烙铁发热器的发热功率为 $x (t)$ （输入），比外界温度高为 $y (t)$ （输出），热导为 $G$ ，热容量为 $C$ ，可列出以下微分方程：
$C\frac{dy(t)}{dt}+Gy(t)=x(t)$
时间常数 $\tau=\frac{C}{G}$
该微分方程的通解 $y(t)=y_0+\frac{1}{C}e^{-t/\tau}\int_{-\infty}^{t}x(t)e^{t/\tau}dt$
先来看看仿真结果

# 参数
C = 20
G = 1
t = np.arange(0, 600, 0.1)
Pf = lambda t: 1 + np.sin(t**2/5000)
P = Pf(t) #发热功率必须为正数

# 生成数据，加入噪声
#pn = P.mean()/50
#P += np.random.normal(0, pn, P.shape)  # 实际功率噪声,odeint无效
# 使用通解进行计算，但数值太大，可能会浮点溢出，不建议使用
T = 1/C*np.exp(-G/C*t)*np.cumsum(P*np.exp(G/C*t))*(t[1]-t[0])
# 使用scipy求解ODE，不能仿真实际功率噪声
#T = odeint(lambda T,t_:(Pf(t_)-T*G)/C, 0, t)[:, 0]
pn = P.mean()/50
P += np.random.normal(0, pn, P.shape)  # 功率测量噪声
tn = T.mean()/50
T += np.random.normal(0, tn, T.shape)  # 温度测量噪声

# 温度测量比功率测量落后(>0)或超前(<0)多少帧
dt = 0
if dt > 0:
    P, T, t = P[:-dt], T[dt:], t[:-dt]
if dt < 0:
    P, T, t = P[-dt:], T[:dt], t[:dt]

# 绘图
plt.plot(t, P, c='r')
plt.twinx()
plt.plot(t, T, c='b')
plt.show()

功率和温度的波形图傅里叶变换一下，看看频谱

Pw = np.fft.rfft(P)
Tw = np.fft.rfft(T)
PwTw = Pw/Tw
w = np.linspace(0, 2*np.pi/(2*(t[1]-t[0])), Pw.shape[0])

n = 30
plt.plot(w[:n], np.log(np.abs(Pw[:n])), c='r')
plt.twinx()
plt.plot(w[:n], np.log(np.abs(Tw[:n])), c='b')
plt.show()

FFT后的频谱图

经过傅里叶变换，变为
$j\omega CY(\omega)+GY(\omega)=X(\omega)$
$(j\omega C+G)Y(\omega)=X(\omega)$
实部 $*Re[Y(\omega)]=Re[X(\omega)]$
虚部 $\omega C*Im[Y(\omega)]=Im[X(\omega)]$
实际测量到的数据是有噪声的
$H^{-1}(\omega)=\frac{X(\omega)}{Y(\omega)}=j\omega C+G$
$G=Re[H^{-1}(\omega)]$
$C=\frac{1}{\omega}Im[H^{-1}(\omega)]$

n = 30
plt.plot(w[:n], PwTw.real[:n], c='r')
plt.twinx()
plt.plot(w[:n], PwTw.imag[:n], c='b')
plt.twiny()
plt.xlim(-1, n+1)  # 数组索引
plt.show()

逆传递函数的实部和虚部
只要X(s)输入啁啾信号，这样可以在实验中测定 $C$ 和 $G$
另外，实际测量到的数据是有噪声的，为了提高分辨率，需要进行最小二乘法
$U)=\sum\left|X_i-(j\omega_iC+G)Y_i\right|^2=\sum\left|\varepsilon_i\right|^2=\sum[Re^2(\varepsilon_i)+Im^2(\varepsilon_i)]$