MATLAB: Matlab code for generating some shapes using signed distance

shapesigned distance

In order to get a good understanding of my question, please, you may first run the matlab code below. This code can generate the shape for rectangle, ellipse and circle using signed distance if you uncomment the portion corresponding to each shape.
Please, I need a code that can give the shapes in the attached picture (Picture_1.jpg) using signed distance. Thank you so much.
%%%Signed distance code for generating shapes
clear all;
close all;
clc;
DomainWidth=2;
DomainHight=1;
ENPC=40;
ENPR=80;
EW = DomainWidth / ENPR; % The width of each finite element.
EH = DomainHight / ENPC; % The hight of each finite element.
M = [ ENPC + 1 , ENPR + 1 ];
[ x,y ] = meshgrid( EW * [ -0.5 : ENPR + 0.5 ] , EH * [ -0.5 : ENPC + 0.5 ]);
[ FENd.x, FENd.y, FirstNdPCol] = MakeNodes(ENPR,ENPC,EW,EH);
LSgrid.x = x(:); LSgrid.y = y(:); % The coordinates of Level Set grid 1
cx = DomainWidth/2;
cy= DomainHight/2;
a = cx;
b = 0.8*cy;
%%Generate circle
tmpPhi= sqrt ( ( LSgrid . x - cx ) .^2 + ( LSgrid . y - cy ) .^2 ) - DomainHight/2;
LSgrid.Phi = -((tmpPhi.')).';
%%Generate ellipse
% tmpPhi= ( (LSgrid . x - cx)/a ) .^2 + (( LSgrid . y - cy)/(0.8*cy) ) .^2 - 1;
% LSgrid.Phi = -((tmpPhi.')).';
%%Generate rectangle
% lower = [cx - 0.5 * DomainWidth, cy - 0.25 * DomainHight];
% upper = [cx + 0.5 * DomainWidth, cy + 0.25 * DomainHight];
% Phi11 =max(LSgrid.x - upper(1), lower(1) - LSgrid . x );
% Phi11 =max(Phi11,LSgrid.y - upper(2));
% Phi11 =max(Phi11,lower(2) - LSgrid . y );
% LSgrid.Phi = -Phi11;
FENd.Phi = griddata( LSgrid.x, LSgrid.y, LSgrid.Phi, FENd.x, FENd.y, 'cubic');
figure(10)
% Figure of the full contour plot of the signed distance of the shape
contourf( reshape( FENd.x, M), reshape(FENd.y , M), reshape(FENd.Phi, M));
axis equal; grid on;drawnow; colorbar;
figure(11)
% Figure of the scaled contour plot showing only the shape
contourf( reshape( FENd.x, M), reshape(FENd.y , M), reshape(FENd.Phi, M), [0 0] );
axis equal; grid on;drawnow; colorbar;
figure(12)
%Figure of the signed distance function
h3=surface(x, y, reshape(-LSgrid.Phi , M + 1)); view([37.5 30]); axis equal; grid on;
set(h3,'FaceLighting','phong','FaceColor','interp', 'AmbientStrength',0.6); light('Position',[0 0 1],'Style','infinite'); colorbar;
%%Code for creating nodes
function [NodesX, NodesY, FirstNdPCol] = MakeNodes(EleNumPerRow,EleNumPerCol,EleWidth,EleHight)
[ x , y ]= meshgrid( EleWidth * [ 0 : EleNumPerRow ], EleHight * [0 : EleNumPerCol]);
FirstNdPCol = find( y(:) == max(y(:)));
NodesX = x(:); NodesY = y(:);
end

Best Answer

I had to read a bit about "signed distance fields" and came across this blogpost
Have your read it? Should not be too hard to implement for your shapes.
I would attempt the following more 'general' algorithm.
  • describe the shape with a set of x-y-coordinates
  • create a mesh
  • find the closest distance to each point in the mesh to the set of x-y-coordinates
If I understand correctly, that is what the "signed distance field" describe, i.e. the closest distance to a shape from any point in the domain.
This is the code for a single horizontal line from [0,0.5] to [1,0.5].
% The "shape"
x=0:0.01:1;
y=ones(size(x)).*.5;
% The "domain"
[X,Y] = meshgrid(0:0.01:1,0:0.01:1);
% Closest point between domain and shape
[X,Y] = meshgrid(0:0.01:1,0:0.01:1);
[id,Z] = dsearchn([x',y'],[X(:),Y(:)]);
scatter3(X(:),Y(:),Z,[],Z)
The only thing you have to change is the sign of those points that are inside of the "shape". This is however fairly simple, using the inpolygon function. Let's try one for the sinusoidal shape:
% Build polygons
x=0:0.01:1;
Y{1} = 0.2.*sin(x/0.05)+0.3
Y{2} = 0.2.*sin(x/0.05)+0.7
py = [Y{1},fliplr(Y{2})]
px = [x,fliplr(x)];
% Define domain
[X,Y] = meshgrid(0:0.01:1,0:0.01:1);
% Calculate distance
[~,Z] = dsearchn([px',py'],[X(:),Y(:)]);
[in]=inpolygon(X(:),Y(:),px,py)
Z(in) = -Z(in)
Z = -Z;
% Plot
surf(X,Y,reshape(Z,size(X)),'edgecolor','none');hold on
p = polyshape(px,py)
plot(p)
axis([0 1 0 1])
.
There seems to be some edge-effect, in particular to the left. Other than that, it looks quite OK.