%Funcao que recebe um Conjunto de pontos de uma resposta impulsiva e devolve o IIR equivalente

function [Z, P, K] = FuncFIR2zpk(y1)
y1 = y1';
%close all;
%clear all;
	% Definir o sistema

	%Num = [0 1]
	%Den = [1 -0.9625]
	%Num = [0 0 1]
	%Den = [1 -1.2115 .2406]	
	%Num = [0 1 -.2 .01]
	%Den = [1 -1.7115 .8464 -.1203]
	%Encontra a resposta impulsiva com n pontos

	%n = 300;
	%y1 = dimpulse(Num, Den, n);

	%load rawdata
	%y1 = Dado0;

%y1=Dado(2,:)';
n = size(y1,1);
%Imprime a resposta
%figure, stem(y1);
%title('Sinal Original');

%Cria o sistema
MatrizA = diag(ones(n-2,1),-1);
MatrizB = zeros(n-1,1); MatrizB(1) = 1;
MatrizC = y1(2:n)';
MatrizD = y1(1);

y2 = dimpulse(MatrizA,MatrizB,MatrizC,MatrizD,1, n);    

%figure, stem(y2);
%title('Sinal Original Convertido para Matriz A a D');

%H = zeros(n-1);
%H(1,1:n-1) = MatrizC(1:n-1);
%H(2,1:n-2) = MatrizC(2:n-1);
%H(3,1:n-3) = MatrizC(3:n-1);
%H(4,1:n-4) = MatrizC(4:n-1);
for indice=1:n-1
   H(indice,1:n-indice) = MatrizC(indice:n-1);
end

[u,s,v] = svd(H);

k=20
Ak = v(2:n-1,1:k)'*v(1:n-2,1:k);
Bk = v(1,1:k)';
Ck = MatrizC*v(1:n-1,1:k);

y3 = dimpulse(Ak,Bk,Ck,0,1, n);    

%figure, stem(y3);
%title('Sinal Recuperado com as Matrizes Ak a Ck (Dk = 0)');

%[Num2, Den2] = ss2tf(Ak, Bk, Ck, 0,1);	 

[Z,P,K] = ss2zp(Ak,Bk,Ck,0, 1);
%y4 = dimpulse(Num2, Den2, n);
%figure, stem(y4);
%title('Sinal Recuperado com a Representacao IIR');

%erro = y4-y1;
%figure, stem(erro);
%title('Sinal de Erro entre o Sinal Original e o Recuperado com a Representacao IIR');

%sing = diag(s);

%figure, stem(sing);
%title('Singular Value');


%ordemP = find(sing>((max(sing)-min(sing))*0.001))