Removing points from meshgrid (mask)

14 views (last 30 days)
António Lúcio
António Lúcio on 20 Nov 2019
Commented: darova on 21 Nov 2019
Hey!
I'm trying to remove points from my mesh grid in order to get a graph similar to this one:
Linear_TOPOG.PNG
Basically I want to remove all the points below the line given in the code below by h1.
clear all; close all; clc % clear variables, close plots, clear screen
%% Constants
D = 60; % Ekman layer depth (m)
rho_0 = 1025; % density (kg/m^3)
f = 0.545/(10^4); % Coriolis parameter
c = (1+i)*pi/D; % c parameter
ug = 0; % zonal geostrophic wind speed
tau_x = 0; % x-component of wind shear stress
tau_y = 0.07; % y-component of wind shear stress
tau = tau_x + i*tau_y; % complex wind shear stress
u0 = (1-i)*pi*tau/(rho_0*f*D); % surface Ekman velocity in the deep ocean limit
hDmax = 2.5; % maximum h/D value
%% Bathymetry
h1 = @(x) -x/1000; % bathymetry profile function
x_grid = linspace(-200000,0,1000); % linearly spaced grid for the horizontal coordinate (offshore distance, x)
z_grid = linspace(-100,0,1000); % linearly spaced grid for the vertical coordinate (depth, z)
h = h1(x_grid); % vector with bathymetry
tmp = find(h/D>hDmax); % reduce offshore grid
xindex = tmp(end);
h = h(xindex:end);
L = length(h);
N = length(z_grid);
zeta = pi*h/D; % vector with zeta function values
%% Horizontally-varying structure functions
alpha = (cosh(zeta).*cos(zeta)).^2 + (sinh(zeta).*sin(zeta)).^2;
S1 = cosh(zeta).*cos(zeta)./alpha;
S2 = sinh(zeta).*sin(zeta)./alpha;
T1 = cosh(zeta).*sinh(zeta)./alpha;
T2 = cos(zeta).*sin(zeta)./alpha;
%% Velocities, transport and stream function
vg = (2*pi*tau_y/(rho_0*f*D)).*(1-S1)./(T1-T2) - ug.*(T1+T2-2.*pi.*h./D)./(T1-T2); % meridional geostrophic wind speed
ugu = h*0 + ug; % adjusting the size of the vector
ugc = complex(ugu,vg); % complex total geostrophic velocity
Uek = abs(tau_y/(rho_0*f)); % Ekman transport
psi = zeros(N,L);
psi = complex(psi,psi);
for jj=1:L
for ii=1:N
psi(ii,jj) = (1/c).*(u0.*(1-(cosh(c.*(z_grid(ii)+h(jj))))./(cosh(c.*h(jj)))) +...
ugc(jj).*sinh(c.*z_grid(ii))./cosh(c.*h(jj)));
end
end
psi = real(psi); % real part of stream function
%% Plotting
X_grid = repmat(-h/D,[N 1]);
Z_grid = repmat(z_grid/D, [L 1]);
Z_grid = Z_grid'
Z = psi/Uek;
contour(X_grid,Z_grid,Z,'ShowText','on','LevelStep',0.1)
The figure produced at the moment is the following:
stream function.jpg
What can I do?

Answers (1)

darova
darova on 20 Nov 2019
Try this simple script
[X,Y] = meshgrid(1:20);
n = [1 -1];
n = n/norm(n); % normal of a line
dist = -7; %
D = n(1)*X + n(2)*Y + dist;
X(D>0) = nan;
Y(D>0) = nan;
surf(X,Y,X*0+1)
The result:
img1.png
  2 Comments
António Lúcio
António Lúcio on 20 Nov 2019
Edited: António Lúcio on 20 Nov 2019
Hey there, thanks for your input, but I was actually thinking a bit more generally than what you sketch in your answer. For instance, for tha graphs below, h is not given by any analytical formula. Just a set of values. How would I eliminate all the values below each h value for each x?
darova
darova on 21 Nov 2019
Use inpolygon to detect points in red region. Then assign NaN
123.png

Sign in to comment.

Categories

Find more on Stress and Strain in Help Center and File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!