{[ promptMessage ]}

Bookmark it

{[ promptMessage ]}

WienHammerParIdent

# WienHammerParIdent - y(k 2 = b(1*y(k 1 b(2*w(k 1 if k>r...

This preview shows page 1. Sign up to view the full content.

%% Parametric identification of wiener-hammerstein system. % Obtain inputs and outputs of model clear all u=normrnd(0,2,1,500); % A white gaussian input sequence u with length %400 0 mean and standard deviation 2 ut=normrnd(0,2,1,200); %input for testing. e=normrnd(0,.2,1,400); % A white gaussian with zero mean and standart de %viation .2 with length 400. it is error term e = zeros(1,400); % this is added after all. actually it should have % been done before N = 400; r = 3; ny = 1; nu = 2; sg = 3; % N: # of training data v = zeros(1,500); w = zeros(1,500); y = zeros(1,500); K = zeros(N-r+1,N-r+1);Ker = zeros(200,200); a = [.8;.5]; b = [.6;.3]; % u = zeros(1,500);u(1,1) = 1; x = zeros(ny+1+nu,500); xr = zeros(ny+1+nu,500); %regression vectors. for k = 1 : 497 v(k+1) = a(1)*v(k) + a(2)*u(k); w(k+1) = sin(v(k+1)); %v(k+1)/(1+v(k+1)); %
This is the end of the preview. Sign up to access the rest of the document.

Unformatted text preview: y(k+2) = b(1)*y(k+1) + b(2)*w(k+1); if( k>r) x(:,k) = [y(k:-1:k-ny) u(k-1:-1:k-nu) ]'; end end figure(1); plot(y,'b'); % [h1,t1] = impz([0 a(2,1)],[1 -a(1,1)]); % vd = conv(h1,u); % wd = sin(vd); %vd./((1+vd)); % % [h2,t2] = impz([0 b(2,1)],[1 -b(1,1)]); % yd = conv(h2,wd); % figure(2) ; plot(yd(1:500),'r'); % hold off % %% construct kernel and other submatrices % Kernel matrix for i = 1:N-r+1 for j = 1:N-r+1 K(i,j) = exp(-((u(i+r-1)-u(j+r-1))^2)); %/(2*sg^2) end end Yf = y(r+2:N+2); Yp = y(r+1:N+1); % construct overall matrix C =800; % Penalty term bigEqMat = [ 0 0 ones(N-r+1,1)' ;... 0 0 Yp ones(N-r+1,1) Yp' b(2)^2*K + (1/C)*eye(N-r+1) ]; rigSide = [0 0 Yf]'; finSolution = bigEqMat\rigSide; d = finSolution(1,1) b1= finSolution(2,1) alph = finSolution(3:end); a...
View Full Document

{[ snackBarMessage ]}

Ask a homework question - tutors are online