%% Series Temporales (ruido coloreado), simulación e identificación
% © 2017, 2022, Antonio Sala Piqueras. Universitat Politècnica de València. Todos los derechos reservados.
% Este código ejecutó en Matlab 2022a
% Presentación en vídeo en: http://personales.upv.es/asala/YT/V/idst.html
% Tabla de Contenidos
% Simulación de Ruido Coloreado
% Identificación de Serie Temporal
% lo pasamos todo a continuo
%
%% Simulación de Ruido Coloreado

s=tf('s');
GruidoColoreado=2/(s+.1);
sgm=norm(GruidoColoreado,2)

%simulemos ruido coloreado con ese ancho de banda:
Ts_sim=0.005;time=0:Ts_sim:1000;
rblanc=randn(size(time))/sqrt(Ts_sim);
[y,x]=lsim(GruidoColoreado,rblanc,time);
sgmy=std(y(10000:end))
plot(time,y)
hold on
plot([0 1000],[-2 2;-2 2]*sgmy,'r'), grid on
hold off

%% Identificación de Serie Temporal

%suponemos que tomamos datos cada 0.05 segundos... el Ts anterior erar para simular un ruido "continuo".
%ademas, datos estacionarios, a partir de 10000 muestras quitar efecto condiciones iniciales.
Ts_ident=10*Ts_sim;
dtaident=iddata(y(10000:10:end),[],Ts_ident);

th=armax(dtaident,[2 2]);
present(th)
th=armax(dtaident,[1 1]);
present(th)
tfnoise=tf(th.c,th.a,Ts_ident)
resid(dtaident,th); grid on
%la excitation tiene desvtip:
sgthid=sqrt(th.NoiseVariance)

%
%% lo pasamos todo a continuo

tfnid=d2c(tfnoise,'tustin')
sigma_ct=sgthid*sqrt(Ts_ident)
estimadoWin_continuo=tfnid*sigma_ct