matlab的M函数PID

%增量式PID控制
clear all; close all;

den=[1 1 0]; num=[1]; Ts=0.5; Td=1; %连续系统对象
sys=tf（numden‘inputdelay‘Td）;
dsys=c2d（sysTs‘zoh‘）;
[dnuma]=tfdata（dsys‘v‘）;
na=length（a）-1; b=dnum（2:na+1）; nb=length（b）-1;
d=Td/Ts+1;

kp=0.4; ki=0.0; kd=1; %PID控制器参数（试凑法）

L=400; %控制步数
uk=zeros（d+nb1）; %输入初值：uk（i）表示u（k-i）
yk=zeros（na1）; %输出初值
ek=zeros（21）; %输出误差初值
yr=10*[ones（L/41）;-ones（L/41）;ones（L/41）;-ones（L/41）]; %期望输出

for k=1:L
    time（k）=k*Ts;
    y（k）=-a（2:na+1）*yk+b*uk（d:d+