基于control的一些python控制处理

1，画波特图

import control as ctr
import matplotlib.pyplot as plt
import numpy as np
import math

def bode(sys,HZ="True",dB="True",color="b",deg="True",w=None):
    ":parameter ：输入一个系统"
    ":returns:如果w没有输入，则输出波特图，而且返回mag和pha以及他们对应地omega"
    "如果有w的输入，则输出频率为w时的pha和mag"
    if w == None:
        mag,pha,omega1 = ctr.bode_plot(sys,Hz=HZ,dB=dB,deg=deg,color=color)
        plt.show()
        return mag,pha,omega1
    if w!= None:
        w1 = [w]
        mag, pha, omega1 = ctr.bode_plot(sys, Hz=HZ, dB=dB, deg=deg, color=color,omega=w1)
        return mag,pha/math.pi * 180

sys = ctr.tf([1],[1,2,1])
mag, pha = bode(sys,w = 5000)
print(mag,pha)

和matlab结果进行比对

2、step函数

def step(sys,time=None):
    '输入：sys，可选项：某一时刻的时间'
    '输出：如果没有时间，则输出图片，以及yout 和 t，这两个列表一一对应'
    if time==None:
        T,yout = ctr.step_response(sys)
        plt.plot(T,yout)
        plt.title = "step response"
        plt.xlabel = "time"
        plt.ylabel = "mag"
        plt.show()
        return yout,T
    if time!=None:
        dt = time/50
        T = np.arange(time-50*dt,time+50*dt,dt)
        T,yout = ctr.step_response(sys,T=T)
        print(T)
        print(yout)
        return yout[49]

sys = ctr.tf([1],[1,2,1])
step(sys)
yout = step(sys,time = 6)
print(yout)

注意，control.step_response函数有个非常坑爹的性质：那就是输出的yout
并不真实反映响应值，而是反映：从0到T，步长为dt的值。如果你的设置是，从0+dt到T+dt，dt为步长，你就会惊喜地发现，这两者输出的yout一模一样。
举个例子：

sys = ctr.tf([1],[1,2,1])
T1 = np.arange(0,7,0.07)
T, yout = ctr.step_response(sys,T1)
print(yout)
T2 = np.arange(1,8,0.07)
T,yout = ctr.step_response(sys,T1)
print(yout)

可以看见一个是从0-7，一个是从1-8，但是两者的响应竟然是一样的

然后是结果对照。由于我没有找到matlab 计算tep具体某个时间点值的函数，所以只能对照一部分结果。

左图是python，右图是matlab

3.impulse函数
这个函数和step一模一样，直接贴代码

def impulse(sys,time=None):
    '输入：sys，可选项：某一时刻的时间'
    '输出：如果没有时间，则输出图片，以及yout 和 t，这两个列表一一对应'
    if time==None:
        T,yout = ctr.impulse_response(sys)
        plt.plot(T,yout)
        plt.title = "impulse response"
        plt.xlabel = "time"
        plt.ylabel = "mag"
        plt.show()
        return yout,T
    if time!=None:
        dt = time/50
        T = np.arange(time-50*dt,time+50*dt,dt)
        T,yout = ctr.impulse_response(sys,T=T)
        print(T)
        print(yout)
        return yout[49]

4.裕量计算

def margin(sys):
    "输入：sys；输出：依次：幅值裕量，相位裕量，幅值=1时的频率，相位等于-180时的频率（注意是频率不是角频率）"
    gain = ctr.stability_margins(sys)[0]
    phase = ctr.stability_margins(sys)[1]
    wc1 = ctr.stability_margins(sys)[3]
    fc = ctr.stability_margins(sys)[4]
    return gain,phase,fc/math.pi/2,wc1/math.pi/2

sys = ctr.tf([20],[0.1,1.1,1,0])
print(margin(sys))
bode(sys)

作为比对

5.根轨迹
简单封装

def rlocus(sys):
    ctr.root_locus(sys)
    plt.show()
    
sys = ctr.tf([20],[0.1,1.1,1,0])
rlocus(sys)

结果显示

6.conv

def conv(a,b):
    "输入了两个列表，输出一个列表"
    alen = len(a)
    blen = len(b)
    c = []
    d = []
    npa = np.array(a)
    npb = np.array(b)
    npc = np.array(c)
    npd = np.array(d)
    for i in range(alen):
        d.append(i*[0] + list(npa[i]*npb)+(alen-i-1)*[0])
    npd = np.array(d)
    npc = np.zeros(int(npd.size/len(a)))
    for i in range(len(a)):
        npc = npc + npd[i]
    return list(npc)
    a = [1.0,2,4]
b = [1,3,0,4]
print(conv(a,b))