main.m

hodge decomposition for 2D vector fields.

%% test_dhhfd.m 
clear all; 
close all;

load px;
load py;

[M,N] = size(px);
%% generate regular rectangular meshes 
%box=2; 
[mv_x,mv_y]=meshgrid(1:M,1:N); 
%% construct a potential function 
zz=-exp(-mv_x.*mv_x-mv_y.*mv_y); 
[uu,vv]=gradient(zz); 
%% synthesize a complex motion field 
mv_u=uu+vv; 
mv_v=vv-uu; 
%% display the synthetic motion field 
figure(1); 
quiver(mv_x,mv_y,px,py);%axis ij; %axis([-box,box,-box,box]); 
%% generate the curl-free and divergence-free potential functions 
[potE,potW]=potential_calculation(mv_x,mv_y,px,py); 
%% display the curl-free potential function E 
figure(2); 
mesh(mv_x,mv_y,potE,'EdgeColor','black'); 
%% display the divergence-free potential function W 
figure(3); 
mesh(mv_x,mv_y,potW,'EdgeColor','black'); 
%% generate the curl-free field 
[uu_E,vv_E]=gradient(potE); 
%% display the curl-free field 
figure(4); 
quiver(mv_x,mv_y,uu_E,vv_E); %axis ij; axis([-box,box,-box,box]);  
%% generate the divergence-free field 
[tmp_u,tmp_v]=gradient(potW); 
uu_W=tmp_v; 
vv_W=-tmp_u 
%% display the divergence-free field 
figure(5); 
quiver(mv_x,mv_y,uu_W,vv_W); %axis ij; axis([-box,box,-box,box]); 
%% calculate the magnitude field of the input motion field 
mag0=sqrt(px.^2+py.^2); 
%% find the largest magnitude of the input motion field 
max_mag0=max(max(mag0)); 
%% add the curl-free field and the divergence-free field 
U1=uu_E+vv_W; 
V1=vv_E-uu_W; 
%% calculate the magnitude field of the added motion field 
mag1=sqrt(U1.^2+V1.^2); 
%% find the largest magnitude of the added motion field 
max_mag1=max(max(mag1)); 
%% calculate the scaling factor, r in the thesis 
ratio=max_mag0/max_mag1; 
%% solve for the harmonic remainder 
ru=mv_u-ratio*U1; 
rv=mv_v-ratio*V1; 
%% display the harmonic remainder 
figure(6); 
quiver(mv_x,mv_y,ru,rv); %axis ij; axis([-box,box,-box,box]);