 # Mesh loft from a pair of polylines

Hi,

Does COMPAS have a mesh loft command?
In my case I have pairs of polylines. In Rhino, I was triangulating polygon using ear clipping algorithm and then creating quad in between them. I am wondering if there is already a ready function.

what is it exactly that you want to do?

Hi,

I would like to input top and bottom polylines and create mesh from them:
Is there is any method that exists on compas, or do I need to write this one myself?

you can use `compas_triangle` for this and if `compas_triangle` is installed in your env you can use it directly via `compas.geometry`

``````from compas.geometry import Polygon
from compas.geometry import constrained_delaunay_triangulation
from compas.datastructures import Mesh, mesh_thicken
from compas_view2.app import App

vertices, faces = constrained_delaunay_triangulation(boundary.points, polygons=[hole.points])

mesh = Mesh.from_vertices_and_faces(vertices, faces)
mesh = mesh_thicken(mesh, 0.1)

viewer = App()
viewer.show()
``````

2 Likes

Thank you very much, delaunay triangulation is what I was looking for.

btw,

fyi, this is now also possible using `compas_cgal`, which has the advantage of being `conda` installable…

Hi,

Getting back on this thing regarding CGAL.

Do you use this sample?
CGAL 5.3 - 2D Triangulations: Triangulation_2/polygon_triangulation.cpp

8.4 Example
CGAL 5.3 - 2D Triangulations: User Manual

Found it, nice implementation:

compas_cgal/triangulations.cpp at main · compas-dev/compas_cgal (github.com)

I really enjoy debugging C++ with COMPAS python, this is much more easier than restarting rhino each time you make a code change.

1 Like

2D constrained delaunay triangulation works without errors (if elements are not intersecting), it is good that I dropped opennurbs and jumped on CGAL and COMPAS wagon. I still need to learn more about package management but it is going well.

Also holes:

1 Like

what do you mean with “if elements are not intersecting”? there is a part of the CGAL triangulation lib that i haven’t wrapped yet that deals with those…

This is the case:

I am wondering if wrapping this to try and catch statement would ignore the crash, at least it works for mesh boolean operations.

It seems I used a bit different code from yours.
This is C++ code.

``````#pragma once
#include "stdafx.h"
#include "compas.h"

namespace CGAL_MeshUtil{

inline void
mark_domains(CDT& ct,
Face_handle start,
int index,
std::list<CDT::Edge>& border) {
if (start->info().nesting_level != -1) {
return;
}
std::list<Face_handle> queue;
queue.push_back(start);
while (!queue.empty()) {
Face_handle fh = queue.front();
queue.pop_front();
if (fh->info().nesting_level == -1) {
fh->info().nesting_level = index;
for (int i = 0; i < 3; i++) {
CDT::Edge e(fh, i);
Face_handle n = fh->neighbor(i);
if (n->info().nesting_level == -1) {
if (ct.is_constrained(e)) border.push_back(e);
else queue.push_back(n);
}
}
}
}
}
//explore set of facets connected with non constrained edges,
//and attribute to each such set a nesting level.
//We start from facets incident to the infinite vertex, with a nesting
//level of 0. Then we recursively consider the non-explored facets incident
//to constrained edges bounding the former set and increase the nesting level by 1.
//Facets in the domain are those with an odd nesting level.
inline void
mark_domains(CDT& cdt) {
for (CDT::Face_handle f : cdt.all_face_handles()) {
f->info().nesting_level = -1;
}
std::list<CDT::Edge> border;
mark_domains(cdt, cdt.infinite_face(), 0, border);
while (!border.empty()) {
CDT::Edge e = border.front();
border.pop_front();
Face_handle n = e.first->neighbor(e.second);
if (n->info().nesting_level == -1) {
mark_domains(cdt, n, e.first->info().nesting_level + 1, border);
}
}
}

inline CGAL::Aff_transformation_3<IK> plane_to_XY_mesh(
IK::Point_3 O0, IK::Plane_3 plane
) {
auto X0 = plane.base1();
auto Y0 = plane.base2();
auto Z0 = plane.orthogonal_vector();
CGAL_VectorUtil::Unitize(X0);
CGAL_VectorUtil::Unitize(Y0);
CGAL_VectorUtil::Unitize(Z0);

// transformation maps P0 to P1, P0+X0 to P1+X1, ...

//Move to origin -> T0 translates point P0 to (0,0,0)
CGAL::Aff_transformation_3<IK> T0(CGAL::TRANSLATION, IK::Vector_3(0 - O0.x(), 0 - O0.y(), 0 - O0.z()));

//Rotate ->
CGAL::Aff_transformation_3<IK> F0(
X0.x(), X0.y(), X0.z(),
Y0.x(), Y0.y(), Y0.z(),
Z0.x(), Z0.y(), Z0.z()
);

return  F0 * T0;

}

inline void mesh_from_polylines(std::vector<CGAL_Polyline>& polylines_with_holes, IK::Plane_3& base_plane, std::vector<int>& top_face_vertices) {
);

//////////////////////////////////////////////////////////////////////////////
//Create Transformation | Orient to 2D
//////////////////////////////////////////////////////////////////////////////
CGAL::Aff_transformation_3<IK> xform_toXY = plane_to_XY_mesh(polylines_with_holes, base_plane);
CGAL::Aff_transformation_3<IK> xform_toXY_Inv = xform_toXY.inverse();

//std::vector<Polygon_2> polygons_2d;
//polygons_2d.reserve((int)(polylines_with_holes.size() * 0.5));

CDT cdt;
for (int i = 0; i < polylines_with_holes.size(); i +=2 ) {

Polygon_2 polygon_2d;
for (int j = 0; j < polylines_with_holes[i].size()-1; j++) {
IK::Point_3 p = polylines_with_holes[i][j].transform(xform_toXY);
polygon_2d.push_back(Point(p.hx(), p.hy()));

}
//polygons_2d.emplace_back(polygon_2d);

//////////////////////////////////////////////////////////////////////////////
//Insert the polygons into a constrained triangulation
//////////////////////////////////////////////////////////////////////////////
cdt.insert_constraint(polygon_2d.vertices_begin(), polygon_2d.vertices_end(), true);
}

//////////////////////////////////////////////////////////////////////////////
//Mark facets that are inside the domain bounded by the polygon
//////////////////////////////////////////////////////////////////////////////
mark_domains(cdt);

int count = 0;
for (Face_handle f : cdt.finite_face_handles()) {
if (f->info().in_domain()) {

++count;

}
}

std::map<CDT::Vertex_handle, int> vertex_index;
int k = 0;
for (auto it = cdt.vertices_begin(); it != cdt.vertices_end(); ++it) {
auto p = it->point();
vertex_index[it] = k;
k++;

}

int number_of_faces = 0;
for (Face_handle f : cdt.finite_face_handles()) {
if (f->info().in_domain()) {
number_of_faces+=3;
}
}

top_face_vertices.reserve(number_of_faces);
for (Face_handle f : cdt.finite_face_handles()) {
if (f->info().in_domain()) {

top_face_vertices.emplace_back(vertex_index[f->vertex(0)]);
top_face_vertices.emplace_back(vertex_index[f->vertex(1)]);
top_face_vertices.emplace_back(vertex_index[f->vertex(2)]);

/*          std::vector<uint32_t> indices;
CGAL::Vertex_around_face_circulator<CGAL::Surface_mesh<CGAL::Exact_predicates_inexact_constructions_kernel::Point_3>> vcirc(out.halfedge(face_index), out), done(vcirc);
do indices.push_back(*vcirc++); while (vcirc != done);
meshRebuilt.SetTriangle(face_index.idx(), indices, indices, indices);*/

}
}

}

}
``````

what i meant with my comment about intersections is that you can tell CGAL what to do in case the input contains intersecting edges, which is a part that i haven’t wrapped yet.

not entirely clear to me how this would be addressed with this snippet. so i don’t really understand what you are asking Got it, do you have any reference to documentation where it shows that?