% Perturbation constained by harmonic shallow-water equations
% Color corresponds to the ampltiude of elevation and contourlines
% are represent the phase.

% domain size
m = 64;
n = 64;

% grid resolution
dx = 2000e3/(m-1);
dy = dx;
[x,y] = ndgrid(dx*[0:m-1],dy*[0:n-1]);


% land-sea mask (true: sea, false: land)
mask = true(m,n);
mask(1:m,[1 n]) = false;
mask(m,1:n) = false;

% bathymetry
h = ones(m,n) * 100;
h(mask == 0) = NaN;

% metric
pm = ones(m,n)/dx;
pn = ones(m,n)/dy;

% parameters
omega = 2*pi/(12*60*60); % angular frequency in 1/s
g = 9.81;
f =  6.8873485167246562E-005;
alpha = 0.0001;

len = 20e3;
Nens = 1;

% number of EOFs retained
k = 50;

% compute contrained perturbations

[Ezeta,Eu,Ev,info] = wce_tides(h,pm,pn,g,f,len,alpha,omega,k,Nens);

%i = sub2ind([m n],2,2);
%[C] = wce_covariance2(info.sv,info.WU,info.lambda,info.WE,info.norm,i);
%[Czeta] = statevector_unpack(info.sv,C);

hold on
contourf(x/1e3,y/1e3,abs(Ezeta),50)
caxis([0 max(abs(Ezeta(:)))]);
shading flat
contour(x/1e3,y/1e3,angle(Ezeta),'k')
xlabel('x (km)');
ylabel('y (km)');
title('Perturbation constained by harmonic shallow-water equations');
hold off
