Skip to content
Snippets Groups Projects
Commit 2fad5299 authored by Felix Kong's avatar Felix Kong
Browse files

kernel_dmd.m predict() tested for heart dataset

parent 8da9ecda
Branches
No related merge requests found
...@@ -132,16 +132,16 @@ classdef kernel_dmd < handle ...@@ -132,16 +132,16 @@ classdef kernel_dmd < handle
X2hatScaled = real(X2hatComplex); X2hatScaled = real(X2hatComplex);
X2hat = obj.ReapplyScaling(X2hatScaled); X2hat = obj.ReapplyScaling(X2hatScaled);
imageSeq = obj.VecSeq2ImSeq(X2hat); imageSeq = uint8(obj.VecSeq2ImSeq(X2hat));
end end
function xhat = predict(obj,initial_condition,num_timesteps) function [xhat,imageSeq] = predict(obj,initial_condition,num_timesteps)
% predicts num_timesteps timesteps from an initial condition % predicts num_timesteps timesteps from an initial condition
% like reconstruct, but more general. % like reconstruct, but more general.
% can be used to predict more timesteps than were in the data % can be used to predict more timesteps than were in the data
tic; tic;
if(isempty(obj.U0) || isempty(obj.S0) || isempty(obj.V0)) if(isempty(obj.Qfull) || isempty(obj.Sfull) || isempty(obj.Vhat))
error('Need to run regression() before running predict().'); error('Need to run regression() before running predict().');
end end
...@@ -149,16 +149,18 @@ classdef kernel_dmd < handle ...@@ -149,16 +149,18 @@ classdef kernel_dmd < handle
xhat = nan(num_timesteps,length(initial_condition)); xhat = nan(num_timesteps,length(initial_condition));
xhat(1,:) = initial_condition; xhat(1,:) = initial_condition;
num_modes = size(obj.Q,1);
for k = 2:num_timesteps for k = 2:num_timesteps
prev_x = xhat(k-1,:); prev_x = xhat(k-1,:);
kernelmat = GaussianKernel(prev_x,obj.X); kernelmat = obj.GaussianKernel(prev_x,obj.X);
kernelmat = kernelmat(1:num_modes)'; % truncate and turn into row vector kernelmat = kernelmat(1:num_modes)'; % truncate and turn into row vector
phi_at_t = kernelmat*obj.Q/obj.S*obj.Vhat; phi_at_t = kernelmat*obj.Q/obj.S*obj.Vhat;
xhat(k,:) = real(phi_at_t*obj.Lambda*obj.Xi); xhat(k,:) = real(phi_at_t*obj.Lambda*obj.Xi);
end end
assert(norm(imag(xhat)) <= 1e-6,'Predicted xhat has non-negligible imaginary part!'); assert(norm(imag(xhat)) <= 1e-6,'Predicted xhat has non-negligible imaginary part!');
xhat = real(xhat); xhat = obj.ReapplyScaling(real(xhat));
imageSeq = uint8(obj.VecSeq2ImSeq(xhat));
predictTime = toc; predictTime = toc;
if(obj.verbose) if(obj.verbose)
......
...@@ -154,6 +154,7 @@ for t = 2:2*size(X,1) ...@@ -154,6 +154,7 @@ for t = 2:2*size(X,1)
xhat(t,:) = real(xnext); xhat(t,:) = real(xnext);
end end
toc; toc;
[xhatDmd,dmd_reconstruction] = D.predict(D.X(1,:),200);
%% Check %% Check
figure(1) figure(1)
...@@ -199,6 +200,8 @@ reconsructionHdl = implay(uint8(reconstruction)); ...@@ -199,6 +200,8 @@ reconsructionHdl = implay(uint8(reconstruction));
set(reconsructionHdl.Parent,'Name',sprintf('Reconstruction: %d modes, sigma2 = %.2e.',r,sigma2)); set(reconsructionHdl.Parent,'Name',sprintf('Reconstruction: %d modes, sigma2 = %.2e.',r,sigma2));
originalHdl = implay(grayImages(:,:,image_start:image_start+num_images)); originalHdl = implay(grayImages(:,:,image_start:image_start+num_images));
set(originalHdl.Parent,'Name','Original data'); set(originalHdl.Parent,'Name','Original data');
dmdReconstructionHdl = implay(dmd_reconstruction);
set(dmdReconstructionHdl.Parent,'Name','DMD');
%% %%
% figure(3); % figure(3);
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment