delete the top and bottom triangular faces of the 3D object

Hi! Is there a way to delete the top and bottom triangular faces of the tube?

nodes = importdata("NODES.mat");
faces = importdata("FACES.mat");

figure
plot3(nodes(:,1),nodes(:,2),nodes(:,3),'k.','Markersize',2);
hold on
trimesh(faces(:,:),nodes(:,1),nodes(:,2),nodes(:,3),'EdgeColor','k','Linewidth',0.1,'Facecolor','b','FaceAlpha',1)
hold off
axis equal

Answer 1

[Fabio Freschi]

By looking at your data I see that the top and bottom triangles have a predominant z direction of the face normals. So I play with this

% your code
nodes = importdata("NODES.mat");
faces = importdata("FACES.mat");
figure
plot3(nodes(:,1),nodes(:,2),nodes(:,3),'k.','Markersize',2);
hold on
trimesh(faces(:,:),nodes(:,1),nodes(:,2),nodes(:,3),'EdgeColor','k','Linewidth',0.1,'Facecolor','b','FaceAlpha',1)
hold off
axis equal

% create triangulation object to easily calculate face normals
TR = triangulation(faces,nodes);
P = incenter(TR);
F = faceNormal(TR);

% plot on the previous plot to see what's happening
hold on  
quiver3(P(:,1),P(:,2),P(:,3), ...
     F(:,1),F(:,2),F(:,3),0.5,'color','r');

% dot product with z directed unit vector
proj = F*[0 0 1].';

% threshold (play with this value)
idx = abs(proj) > 0.6;

% remove faces
faces(idx,:) = [];

% new plot without faces
figure
trimesh(faces(:,:),nodes(:,1),nodes(:,2),nodes(:,3),'EdgeColor','k','Linewidth',0.1,'Facecolor','b','FaceAlpha',1)

Hope it helps

[Alberto Acri]

Thank you for your answer @Fabio Freschi. I have several geometries and your code doesn't always work. Is it possible to apply your code in restricted areas?

For example, for this new geometry, I would like to be able to delete only the triangles that lie on the magenta planes (or possibly in the area between the green/red planes).

I thank you if you can help me!

nodes = importdata("nodes_geometry.txt");
faces = importdata("faces_geometry.txt");

P1 = [-24.4053 -9.7273 39.6897];
P1_sup = [-24.24285 -9.75115 38.98785];
P1_inf = [-24.56465 -9.73375 40.36425];
V1 = importdata("V1_geometry.mat");
V1_sup = importdata("V1_geometry_SUP.mat");
V1_inf = importdata("V1_geometry_INF.mat");

P2 = [-25.0631 -13.8062 50.2981];
P2_sup = [-24.8931 -14.06965 50.88635];
P2_inf = [-25.19075 -13.4698 49.7241];
V2 = importdata("V2_geometry.mat");
V2_sup = importdata("V2_geometry_SUP.mat");
V2_inf = importdata("V2_geometry_INF.mat");

figure
plot3(nodes(:,1),nodes(:,2),nodes(:,3),'k.','Markersize',2);
hold on
trimesh(faces(:,:),nodes(:,1),nodes(:,2),nodes(:,3),'EdgeColor','k','Linewidth',0.1,'Facecolor','b','FaceAlpha',0.5)

plot3(P1(:,1),P1(:,2),P1(:,3),'m.','Markersize',20);
patch(V1(:,1), V1(:,2),V1(:,3),'m','FaceAlpha',1);

plot3(P1_sup(:,1),P1_sup(:,2),P1_sup(:,3),'g.','Markersize',10);
patch(V1_sup(:,1), V1_sup(:,2),V1_sup(:,3),'g','FaceAlpha',0.5);

plot3(P1_inf(:,1),P1_inf(:,2),P1_inf(:,3),'r.','Markersize',10);
patch(V1_inf(:,1), V1_inf(:,2),V1_inf(:,3),'r','FaceAlpha',0.5);


plot3(P2(:,1),P2(:,2),P2(:,3),'m.','Markersize',20);
patch(V2(:,1), V2(:,2),V2(:,3),'m','FaceAlpha',1);

plot3(P2_sup(:,1),P2_sup(:,2),P2_sup(:,3),'g.','Markersize',10);
patch(V2_sup(:,1), V2_sup(:,2),V2_sup(:,3),'g','FaceAlpha',0.5);

plot3(P2_inf(:,1),P2_inf(:,2),P2_inf(:,3),'r.','Markersize',10);
patch(V2_inf(:,1), V2_inf(:,2),V2_inf(:,3),'r','FaceAlpha',0.5);


hold off
axis equal

[Fabio Freschi]

I don't have nodes_geometry.txt and faces_geometry.txt to try.

I think about many possibilities

1. if you have the magenta plane, you know the normal to the elements to delete. Use the previous procedure to select those elements and delete it. In this case the result is less sensible to the threshold.
2. delete the elements whose barycenter lays on the magenta plane (with tolerance)
3. You can combine the previous methods using a logic and
4. identify the border of your mesh (check featureEdge), Provide a point on the surface. Find the triangle containing that point (with tolerance). Select the triangles connected to the strating one and repeat this operation until you reach one triangle that is on the border (i.e. that has an edge identified by featureEdge)

[Alberto Acri]

Thanks for your answer @Fabio Freschi! I have attached the files here!

Could you tell me how to proceed to eliminate the triangles lying on the magenta plane?

[Fabio Freschi]

You can use the method I previpusly suggested, using the info regarding the unit vector orthogonal to the magenta planes. Maybe you can do this in a smarter way, starting from the knowledge of the uquations of magenta planes.

clear variables, close all

nodes = importdata("nodes_geometry_1.mat");
faces = importdata("faces_geometry_1.mat");

P1 = [-24.4053 -9.7273 39.6897];
V1 = importdata("V1_geometry.mat");

P2 = [-25.0631 -13.8062 50.2981];
V2 = importdata("V2_geometry.mat");

figure, axis equal, hold on
plot3(nodes(:,1),nodes(:,2),nodes(:,3),'k.','Markersize',2);
trimesh(faces(:,:),nodes(:,1),nodes(:,2),nodes(:,3),'EdgeColor','k','Linewidth',0.1,'Facecolor','b','FaceAlpha',0.5)

plot3(P1(:,1),P1(:,2),P1(:,3),'m.','Markersize',20);
patch(V1(:,1), V1(:,2),V1(:,3),'m','FaceAlpha',1);

plot3(P2(:,1),P2(:,2),P2(:,3),'m.','Markersize',20);
patch(V2(:,1), V2(:,2),V2(:,3),'m','FaceAlpha',1);


% create triangulation object to easily calculate face normals
TR = triangulation(faces,nodes);
P = incenter(TR);
F = faceNormal(TR);

% get unit vector from V1
u1 = cross(V1(1,:)-P1,V1(2,:)-P1);
u1 = u1/sqrt(sum(u1.^2,2));
quiver3(P1(1),P1(2),P1(3),u1(1),u1(2),u1(3));

% get unit vector from V2
u2 = cross(V2(1,:)-P2,V2(2,:)-P2);
u2 = u2/sqrt(sum(u2.^2,2));
quiver3(P2(1),P2(2),P2(3),u2(1),u2(2),u2(3));

% find normals oriented like u1 (dot product)
idx1 = abs(F*u1.') > 0.9999;

% find normals oriented like u2 (dot product)
idx2 = abs(F*u2.') > 0.9999;

% remove bottom lid
faces(idx1 | idx2,:) = [];

% new plot without faces
figure, axis equal
trimesh(faces(:,:),nodes(:,1),nodes(:,2),nodes(:,3),'EdgeColor','k','Linewidth',0.1,'Facecolor','b','FaceAlpha',1)

[Alberto Acri]

Thanks! A doubt that has come to me now. Does the code remove the faces on the magenta plane + the triangular faces that are parallel to that plane?

[Fabio Freschi]

yes. This is why I suggested you also to check barycenters laying on the magenta plane (item 2) or combine the two methods (item 3). To check points on the magenta plane, you can use the eqaution of the plane and check if the barycenters fit this equation. In the previous code, the index of faces to remove are

% find barycenters laying on the magenta plane
idx1 = abs(sum(P.*u1,2)-sum(u1.*P1,2)) < 0.001;

% find barycenters laying on the magenta plane
idx2 = abs(sum(P.*u2,2)-sum(u2.*P2,2)) < 0.001;

Note that also in this case, if you have a complex mesh and another barycenter lays on the plane "far" from your region of interst, that face will be removed

[Alberto Acri]

Thank you. This method seems to work better than the first method.

But to combine the two methods, do I need to combine these two?

idx1_normals = abs(F*u1.') > 0.9999;
idx1_barycenters = abs(sum(P.*u1,2)-sum(u1.*P1,2)) < 0.001;
idx1 = idx1_normals + idx1_barycenters;  % ??

idx2_normals = abs(F*u2.') > 0.9999;
idx2_barycenters = abs(sum(P.*u2,2)-sum(u2.*P2,2)) < 0.001;
idx2 = idx2_normals + idx2_barycenters;  % ??

thanks in advance!

[Fabio Freschi]

the idxs are logical vectors, so they must be combined using logical operator, in this case &

[Alberto Acri]

Is this correct in this way?

idx1 = idx1_normals & idx1_barycenters;
idx2 = idx2_normals & idx2_barycenters;

[Fabio Freschi]

give it a go and see what happens!

[Alberto Acri]

I think it doesn't work :/

[Fabio Freschi]

% find normals oriented like u1 (dot product)
idx1 = abs(F*u1.') > 0.9999;

% find normals oriented like u2 (dot product)
idx2 = abs(F*u2.') > 0.9999;

% find barycenters laying on the magenta plane
ipln1 = abs(sum(P.*u1,2)-sum(u1.*P1,2)) < 0.001;

% find barycenters laying on the magenta plane
ipln2 = abs(sum(P.*u2,2)-sum(u2.*P2,2)) < 0.001;

% remove bottom lid
faces((idx1 & ipln1) | (idx2 & ipln2),:) = [];

[Alberto Acri]

Thank you!