%図14.11, 図14.12, 図14.13
%
clear all; close all;

%図14.11のプロット

%制御対象のパラメータを与える
K = 1; %K=1
omegan = 0.1; %ω=0.1
zeta = 0.2; %ζ=0.2

%制御対象の分子・分母多項式を与える
num = [ 0 0 K * omegan^2 ]; %制御対象の分子多項式
den = [ 1 2*zeta*omegan omegan^2 ]; %制御対象の分母多項式

%コントローラC_0の分子・分母多項式を与える
numc0 = [ 1 ]; %C_0の分子多項式
denc0 = [ 1 ]; %C_0の分母多項式

%コントローラC_1=K_pの分子・分母多項式を与える
numc1 = [ 100 ]; %C_1の分子多項式
denc1 = [ 1 ]; %C_1の分母多項式

%位相遅れコントローラの分子・分母多項式を与える
omega1 = 0.1; %位相遅れコントローラのパラメータ
numPLG = [ 1 omega1 ]; %位相遅れコントローラの分子多項式
denPLG = [ 1 0 ]; %位相遅れコントローラの分母多項式

%位相進みコントローラのパラメータと分子・分母多項式を与える
omega3 = 8; %位相進みコントローラのパラメータ
omega4 = 0.5; %位相進みコントローラのパラメータ
numPLE = [ omega3 omega3 * omega4 ]; %位相進みコントローラの分子多項式
denPLE = [ omega4 omega4 * omega3 ]; %位相進みコントローラの分母多項式

%Roll offコントローラのパラメータと分子・分母多項式を与える
omega5 = 10; %ω_5=10
zeta5 = 0.5; %ζ_5=0.5
num5 = [omega5^2]; %Roll offコントローラの分子多項式
den5 = [1, 2 * zeta5 * omega5, omega5^2]; %Roll offコントローラの分母多項式

%制御対象とコントローラの伝達関数表現を与える
sys = tf( num, den ); %制御対象の伝達関数表現
c0 = tf( numc0, denc0 ); %C_0の伝達関数表現
c1 = tf( numc1, denc1 ); %C_1の伝達関数表現
cPLG = tf( numPLG, denPLG ); %位相遅れコントローラの伝達関数表現
cPLE = tf( numPLE, denPLE ); %位相進みコントローラの伝達関数表現
cont_roll = tf(num5, den5); %Roll offコントローラの伝達関数表現

%開ループ伝達関数を求める
sysL3 = c1*cPLG*cPLE*sys; %L_3(s) = P(s)C_3(s)の場合の開ループ伝達関数
sysL4 = c1*cPLG*cPLE*cont_roll*sys; %L_4(s) = P(s)C_4(s)の場合の開ループ伝達関数

%制御対象とコントローラのフィードバック結合を求める
sysc3 = feedback( sysL3, 1 ); %L_3(s) = P(s)C_3(s)の場合のフィードバック結合
sysc4 = feedback( sysL4, 1 ); %L_4(s) = P(s)C_4(s)の場合のフィードバック結合

%図14.16のプロット

figure(1) %図のウィンドウを開く

%時間変数の定義
t = 0:0.01:10; %0から10まで0.01刻み

%ステップ信号の定義
us = ones( size( t ) );

%観測ノイズの定義
ns = 0.3 * sin( 25*t);

%ステップ応答の計算
[ y3, t3 ] = lsim( sysc3, us + ns, t ); %L_3(s) = P(s)C_3(s)の場合のステップ応答（lsimによる計算）
[ y4, t4 ] = lsim( sysc4, us + ns, t ); %L_4(s) = P(s)C_4(s)の場合のステップ応答（lsimによる計算）

%図14.16(a)のプロット
subplot(1,2,1); %複数の図を並べるためのコマンド．1行2列の1列目という意味
plot(t3,y3); %ステップ応答をプロット
xlim([0,10]) %横軸（時間軸）の範囲の指定
ylim([0.0 1.5]) %縦軸の範囲の指定
xticks([0 2 4 6 8 10]) %横軸の目盛りの値の設定
yticks([0.0 0.5 1.0 1.5]) %縦軸の目盛りの値の設定
grid; %罫線を表示
xlabel('time  t [s]'); %横軸のラベル表示
ylabel('y(t)') %縦軸のラベル表示
title('C_3(s)'); %タイトルの表示

%図14.16(b)のプロット
subplot(1,2,2); %複数の図を並べるためのコマンド．1行2列の2列目という意味
plot(t4,y4); %ステップ応答をプロット
xlim([0,10]) %横軸（時間軸）の範囲の指定
ylim([0.0 1.5]) %縦軸の範囲の指定
xticks([0 2 4 6 8 10]) %横軸の目盛りの値の設定
yticks([0.0 0.5 1.0 1.5]) %縦軸の目盛りの値の設定
grid; %罫線を表示
xlabel('time  t [s]'); %横軸のラベル表示
ylabel('y(t)') %縦軸のラベル表示
title('C_4(s)') %タイトルの表示

%図14.17

%コントローラC(4)を与える
c4 = c1*cPLG*cPLE*cont_roll;

%角周波数の範囲を指定
w = logspace(-2, 2, 1000); %対数的に等間隔なベクトルの生成(10^{-2}から10^{2}で1000点)

%図14.17のプロット

figure(2) %図のウィンドウを開く

%図14.17(a)のプロット
[ gainc4 phasec4 wc4 ] = bode( c4, w ); %ゲインと位相
gainc4 = gainc4(:); %C_4(s)のゲインの配列を1次元ベクトルに変換
gainc4dB = 20*log10(gainc4); %C_4(s)のゲインをデシベル値に変換
phasec4 = phasec4(:); %C_4(s)の位相の配列を1次元ベクトルに変換

subplot(2,2,1) %複数の図を並べるためのコマンド．2行1列の1行1列目という意味
semilogx(wc4, gainc4dB); %ゲイン線図をプロット
xlim([10^(-2),10^(2)]) %横軸（角周波数）の範囲の指定
ylim([0 70]) %縦軸の範囲の指定
xticks([10^(-1) 10^(0) 10^(1)]) %横軸の目盛りの値の設定
yticks([0 20 40 60]) %縦軸の目盛りの値の設定
grid; %罫線を表示
xlabel('Frequency (rad/s)'); %横軸のラベル表示
ylabel('Gain [dB]'); %縦軸のラベル表示
title('C_4(s)') %タイトルの表示

subplot(2,2,3) %複数の図を並べるためのコマンド．2行1列の2行1列目という意味
semilogx(wc4, phasec4); %位相線図をプロット
xlim([10^(-2),10^(2)]) %横軸（角周波数）の範囲の指定
ylim([-90 50]) %縦軸の範囲の指定
xticks([10^(-1) 10^(0) 10^(1)]) %横軸の目盛りの値の設定
yticks([-90 -45 0 45]) %縦軸の目盛りの値の設定
grid; %罫線を表示
xlabel('Frequency (rad/s)'); %横軸のラベル表示
ylabel('Phase [deg]'); %縦軸のラベル表示

%図14.17(b)のプロット
[ gainL4 phaseL4 wL4 ] = bode( sysL4, w ); %ゲインと位相
gainL4 = gainL4(:); %L_4(s) = P(s)C_4(s)の場合のゲインの配列を1次元ベクトルに変換
gainL4dB = 20*log10(gainL4); %L_4(s) = P(s)C_4(s)の場合のゲインをデシベル値に変換
phaseL4 = phaseL4(:); %L_4(s) = P(s)C_4(s)の場合の位相の配列を1次元ベクトルに変換

subplot(2,2,2) %4複数の図を並べるためのコマンド．2行1列の1行2列目という意味
semilogx(wL4, gainL4dB); %ゲイン線図をプロット
xlim([10^(-2),10^(2)]) %横軸（角周波数）の範囲の指定
ylim([-60 60]) %縦軸の範囲の指定
xticks([10^(-1) 10^(0) 10^(1)]) %横軸の目盛りの値の設定
yticks([-40 0 40]) %縦軸の目盛りの値の設定
grid; %罫線を表示
xlabel('Frequency (rad/s)'); %横軸のラベル表示
ylabel('Gain [dB]'); %縦軸のラベル表示
title('L_4(s) = P(s)C_4(s)') %タイトルの表示

subplot(2,2,4) %複数の図を並べるためのコマンド．2行1列の2行2列目という意味
semilogx(wL4, phaseL4); %位相線図をプロット
xlim([10^(-2),10^(2)]) %横軸（角周波数）の範囲の指定
ylim([-225 0]) %縦軸の範囲の指定
xticks([10^(-1) 10^(0) 10^(1)]) %横軸の目盛りの値の設定
yticks([-180 -90 0]) %縦軸の目盛りの値の設定
grid; %罫線を表示
xlabel('Frequency (rad/s)'); %横軸のラベル表示
ylabel('Phase [deg]'); %縦軸のラベル表示

Published with MATLAB® R2016b