clear all 
close all 
clc

%%% Capitulo I 
%%%% Solucao de equacoes diferenciais - Pag 38

% 
s = tf('s');
H0 = (s+3)/(s^2+3*s+2);
H0 = zpk(H0); % Coloca no formato de ganho, polos e zeros
H = (s-3)/(s^2+3*s+2);

t = 0:0.01:10;
h0 = impulse(H0,t);

% Para g(t) degrau unitário
hcg = step(H,t); % h conv g

y = h0+hcg;
plot(t,y)

%hat y(s)  = H0(s)+H(s)g(s) = Ht(s) % 
Ht = (s^2+4*s-3)/(s^3+3*s^2+2*s);
yt = impulse(Ht,t); 
hold on
plot(t,yt)

% Para g(t) = sin(2*t);
g = sin(2*t);
hcg = lsim(H,g,t);
y = h0+hcg;
figure
plot(t,y)

% Para g(t) rampa unitária
hcg = step(H/s,t);
y = h0+hcg;
figure
plot(t,y)
%%%%%%%%%%%%%%%
%%% Tópico: Representacao em Espaco de EStado
% Exemplo dado em aula
% Sistema na forma canonica controlavel
A = [0 1 0;0 0 1;-6 -11 -7];
B = [0;0;1];
C = [16 8 1];
D = 0;

% Obtencao da funcao de transferencia
H = tf(ss(A,B,C,D));

% O mesmo resultado pode ser obtido usando os comandos 
[Num, Den] = ss2tf(A,B,C,D);
H = tf(Num,Den);

% Comando para obter a realizacao minima, ja considerando cancelamentos
% de polos e zeros
H = minreal(H);


% Obtencao das matrizes (A,B,C,D) a partir da funcao de transferencia
Num = [1 8 16];
Den = [1 7 11 6];

[A,B,C,D] = tf2ss(Num,Den);
%%%%%%%%%%%%%%%
% Reposta em frequencia 
% Exemplo pag 52 
clear all 
close all 
clc 
s = tf('s');
H = (6*s+100)/(s^2+6*s+100);
t = 0:0.01:4;
g = (1/4)*sin(10*t);

% Obtencao da saida y(t)
y = lsim(H,g,t);

% Calculo da saida em regime permanente yperm(t)
a = 10*sqrt(-1);
H = (6*a+100)/(a^2+6*a+100);
mH = abs(H);
argH = phase(H);
yperm = mH*(1/4)*sin(10*t+argH);
plot(t,y,t,yperm)
hold on 
plot(t,g)

%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Motor de Corrente Continua
%%%%%%%%%%%%%
% Modelagem e resposta dinamica do motor CC
% pag 92
clear all 
close all 
clc
c=2;K=1/2;Jc=10;Jm=1;b=1.0;L=1;R=1;
% Velocidade desejada
vperm = 50;
JT = Jc+Jm*(c^2);

% Calculo da tensao necessaria para obter a velocidade desejada
V0 = vperm*(R*b+c^2*K^2)/(c*K)

s = tf('s')
G = c*K/((JT*s+b)*(L*s+R)+c^2*K^2)

b = 0.5;
G1 = c*K/((JT*s+b)*(L*s+R)+c^2*K^2);

b = 1.5;
G2 = c*K/((JT*s+b)*(L*s+R)+c^2*K^2);
t = 0:0.01:40;
y = step(V0*G,t);
y1 = step(V0*G1,t);
y2 = step(V0*G2,t);

figure('Position', [10 10 750 270])
plot(t,y,'b','LineWidth',2)
hold on 
plot(t,y1,'--r','LineWidth',2)
hold on 
plot(t,y2,':g','LineWidth',2)
xlabel('$t$','FontSize',20,'Interpreter','latex');
ylabel('$\nu(t)$','FontSize',20,'Interpreter','latex');










