Mesh loft from a pair of polylines


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?


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 import App

boundary = Polygon.from_sides_and_radius_xy(4, 2)
hole = Polygon.from_sides_and_radius_xy(3, 1)

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()


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


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


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 (

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) {
        std::list<Face_handle> queue;
        while (!queue.empty()) {
            Face_handle fh = queue.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();
            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();

        // 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[0][0], 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()));


            //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
        int count = 0;
        for (Face_handle f : cdt.finite_face_handles()) {
            if (f->info().in_domain()) {


        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;


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

        for (Face_handle f : cdt.finite_face_handles()) {
            if (f->info().in_domain()) {

      /*          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[0], indices[1], indices[2]);*/




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 :slight_smile:

Got it, do you have any reference to documentation where it shows that?