Refine delaunay triangulation mesh

I have the code wich generates this delaunay triangulation. I need to add more triangles inside this object to get it as uniform as possible. I wanted to refine my mesh in order to remove those chunky triangles inside and get a larger number of triangles.

Note that my object has 3 holes inside it, and obviously, I dont want to fill them with triangles

How do I get a refined delaunay triangulation mesh? Without using the alphashape command

Answer 1

[Kevin Holly]

I was able to increase the density using the subtri() function upload by Kevin Moerman to the File Exchange.

ImgClosed = imread('test.png');
ImgClosed=im2bw(ImgClosed);
[B,L,N,A] = bwboundaries(ImgClosed); 
figure
imshow(ImgClosed)
hold on; 
for k = 1:N
    % Boundary k is the parent of a hole if the k-th column of the adjacency matrix A contains a non-zero element
    if (nnz(A(:,k)) > 0)
        boundary = B{k};
        redline = plot(boundary(:,2),boundary(:,1),'r','LineWidth',2);
        % Loop through the children of boundary k
        for l = find(A(:,k))'
            boundary = B{l};
            greenline(l,:) = plot(boundary(:,2),boundary(:,1),'g','LineWidth',2);
        end
    end
end

figure
axes
plot(redline.XData,redline.YData,'r','LineWidth',2)
hold on
x=redline.XData;
y=redline.YData;
c = [(1:length(x)-1).' (2:length(x)).'];
for jj = 1:size(greenline,1)
    try
        plot(greenline(jj).XData,greenline(jj).YData,'g','LineWidth',2);
        x_start = length(x);
        x=[x greenline(jj).XData];
        y=[y greenline(jj).YData];
        c = [c;(x_start+(1:length(greenline(jj).XData)-1)).' (x_start+(2:length(greenline(jj).XData))).'];
    catch
    end
end
DT = delaunayTriangulation(x',y',c);

Warning: Duplicate data points have been detected and removed.
The Triangulation indices and constraints are defined with respect to the unique set of points in delaunayTriangulation.

triplot(DT)

figure
DT = delaunayTriangulation(x',y',c);

Warning: Duplicate data points have been detected and removed.
The Triangulation indices and constraints are defined with respect to the unique set of points in delaunayTriangulation.

IO = isInterior(DT);
triplot(DT(IO, :),DT.Points(:,1), DT.Points(:,2),'LineWidth', 2)

h=trimesh(DT(IO, :),DT.Points(:,1), DT.Points(:,2),zeros(size(DT.Points(:,2))),'LineWidth', 2);

[Fs,Vs]=subtri(h.Faces,h.Vertices,2);% Note, you can increase the density by increase the value of the last input
figure
subplot(1,2,1);patch('Faces',h.Faces,'Vertices',h.Vertices,'FaceAlpha',0.5,'EdgeColor','k','LineWidth',2);
hold on;
subplot(1,2,2);patch('Faces',Fs,'Vertices',Vs,'FaceAlpha',0.5,'EdgeColor','k','LineWidth',2);

[Kevin Moerman]

@Kevin Holly @Bruno Neto thanks for checking out subtri, note the latest version of it is part of my GIBBON project (see also: https://www.mathworks.com/matlabcentral/fileexchange/48208-gibboncode-gibbon and https://github.com/gibbonCode/GIBBON/blob/master/lib/subtri.m).

@Bruno Neto if you want it as uniform as possible subtri is not great in this case. It would be better to use: https://www.gibboncode.org/html/HELP_regionTriMesh2D.html. You'd first need the boundary edges using https://www.gibboncode.org/html/HELP_patchBoundary.html then group each to get your inner edges and outer edges e.g. https://www.gibboncode.org/html/HELP_tesgroup.html then convert each edge list to a curve using https://www.gibboncode.org/html/HELP_edgeListToCurve.html .

If you use GIBBON for research, please cite it using: https://doi.org/10.21105/joss.00506

Let me know if you have questions: [e-mail]