function [A,B,D,E,x0]=mydsr(Y,U,L,g,J);
%  [A,B,D,E,x0]=mydsr(Y,U,L,g,J);

[Ny,ny]=size(Y); 
[Nu,nu]=size(U);
N=min(Ny,Nu);

% Genererer datamatriser
[Uf,Yf,Up,Yp,YL,UL]=mat4oppg13(Y,U,L,g,J);

Yt = [Uf;Up;Yp];
Zd = Yf*Yt'*pinv(Yt*Yt')*Yt;                 % Zd_(J/L+1)
Z  = Zd-Zd*Uf'*pinv(Uf*Uf')*Uf;              % Z_(J/L+1)

% Oppgave 13d);
Z0 = Z(1:L*ny,:);
[U,S,V]=svd(Z0);
diag(S)'
n=1;
n=dread('Specify system order: n=',n);

U1=U(:,1:n); S1=S(1:n,1:n); V1=V(:,1:n);
O=U1
D=O(1:ny,:);

% Oppgave 13e)
Z1=Z(ny+1:(L+1)*ny,:);
A=O'*Z1*V1*pinv(S1)

% Oppgave 13f)
Atilde = O*A*pinv(O'*O)*O'

% Oppgave 13g)
Zd0=Zd(1:L*ny,:);
Zd1=Zd(ny+1:(L+1)*ny,:);
Btilde = (Zd1-Atilde*Zd0)*Uf'*pinv(Uf*Uf')

% Oppgave 13h) 
% Beregning av E.
E=Btilde((L-1)*nu+1:L*nu,(L+g-1)*nu+1:(L+g)*nu)

% Beregning av B.
Ei=Btilde(:,(L+g-1)*nu+1:(L+g)*nu);                        % Siste colonne i Btilde.
Hd=zeros(L*ny,(L+g-1)*nu); 
for ic=1:L+g-1
    Hd(:,(L+g-1-ic)*nu+1:(L+g-ic)*nu)=Ei;
    Bc=Btilde(:,(L+g-1-ic)*nu+1:(L+g-ic)*nu);               % Colonne i i Btilde
    Ei=Bc+Atilde*Ei;
end
B=pinv(O'*O)*O'*Ei

% Initialverdi, x0.
x0=pinv(O'*O)*O'*(YL(1:L*ny,1)-Hd*UL(1:L*nu,1));