function    [r0k,omega,Tk,G,Gr]=DCmotor(J,b,K,R,L)

%Vstupni hodnoty:
J=0.02;    %moment setrvacnosti
b=0.2;     %koeficient treni
K=0.2;     %konstanta umernosti
R=2;       %odpor
L=0.4;     %indukcnost


%Prvni vypocet a graf
num=K   %citatel tf
den=[J*L,J*R+L*b,b*R+K^2 0]  %jmenovatel tf

motor=tf(num,den)  %tf-prenosova funkce
motor=minreal(motor)  %rozklad jmenovatele

zpk(motor)  %prenosova funkce ve forme polu a nul

figure, rlocus(motor), grid   %graf polu a nul
axis([-13 2 -9 9]);


%Druhy vypocet a graf:
r0k=(J*R+b*L)*(b*R+K^2)/(J*L*K)     %kriticke zesileni
omega=sqrt((K*r0k)/(J*R+b*L))      %uhlova rychlost
Tk=2*pi/omega                      %perioda


%Optimalni nastaveni pro danou r0k a Tk pro P a PID regulator
%Vypocet a vykresleni prubehu pro P
Kp=0.5*r0k
G=[Kp]
P=feedback(G*motor,1)  %zpetna vazba
%Metoda serizeni Ziegler-Nichols pro PID
Kp=0.6*r0k
Ki=2*Kp/Tk
Kd=Kp/8*Tk
Gr=tf([Kd Kp Ki],[1 0]) %prenosova funkce
PID=feedback(Gr*motor,1)  %zpetna vazba
figure, step(P,PID)  %vykresleni krivek pro P a PID do grafu


end