図6.1のプロット

#計算とグラフプロットに必要なモジュールの読み込み
import numpy as np
from control import matlab
from matplotlib import pyplot as plt
from scipy import arange 

#2次遅れ系のパラメータを与える
K = 1.0 #K=1
omegan = 1.0 #ω=1
zeta1 = 0.1 #ζ=0.1
zeta2 = 1.0 #ζ=1
zeta3 = 2.0 #ζ=2

#伝達関数の分子・分母多項式を与える
num = [0, 0, K * omegan**2] #分子多項式
den1 = [1, 2 * zeta1 * omegan, omegan**2] #ζ=0.1の場合の分母多項式
den2 = [1, 2 * zeta2 * omegan, omegan**2] #ζ=1の場合の分母多項式
den3 = [1, 2 * zeta3 * omegan, omegan**2] #ζ=2の場合の分母多項式

#伝達関数表現を与える
sys1 = matlab.tf(num, den1) #ζ=0.1の場合の伝達関数表現
sys2 = matlab.tf(num, den2) #ζ=1の場合の伝達関数表現
sys3 = matlab.tf(num, den3) #ζ=2の場合の伝達関数表現

#時間変数の定義
t = arange(0, 20, 0.01) #0から20まで0.01刻み

#インパルス応答の計算
y1, t1 = matlab.impulse(sys1, t) #ζ=2の場合のインパルス応答
y2, t2 = matlab.impulse(sys2, t) #ζ=2の場合のインパルス応答
y3, t3 = matlab.impulse(sys3, t) #ζ=2の場合のインパルス応答

#図6．1のプロット
plt.plot(t1, y1, label = "zeta = 0.1") #インパルス応答をプロット
plt.plot(t2, y2, label = "zeta = 1") #インパルス応答をプロット
plt.plot(t3, y3, label = "zeta = 2") #インパルス応答をプロット
plt.xlim([0,20]) #横軸（時間軸）の範囲の指定
plt.xticks([0,5,10,15,20]) #横軸の表示の指定
plt.yticks([-0.8,-0.6,-0.4,-0.2,0,0.2,0.4,0.6,0.8,1.0]) #縦軸の表示の指定
plt.grid(color='gray') #罫線を灰色で表示
plt.xlabel("time   t[s]") #横軸のラベル表示
plt.ylabel("y(t)") #縦軸のラベル表示
plt.legend() #凡例の表示
plt.show() #グラフの表示