VTK Learning Thirty - vtkPolyData to CGAL Surface_mesh

VTK Learning Thirty - vtkPolyData to CGAL Surface_mesh

Description

将vtkPolyData转换成CGAL::Surface_mesh,并利用CGAL的布尔运算对数据进行处理,
然后再将处理后的CGAL::Surface_mesh转换成vtkPolyData进行可视化渲染。

定义 CGAL::Surface_mesh

typedef CGAL::Exact_predicates_inexact_constructions_kernel EPICK;
typedef EPICK::Point_3 Point_3;
typedef CGAL::Surface_mesh<Point_3> SMesh;
typedef boost::property_traits<boost::property_map<SMesh,CGAL::vertex_point_t>::type>::value_type Point;

vtkPolyData to CGAL::Surface_mesh

bool vtkPointSet_to_polygon_mesh(vtkPointSet *poly_data, SMesh &tmesh)
{
    typedef typename boost::property_map<SMesh, CGAL::vertex_point_t>::type VPMap;
    typedef typename boost::property_map_value<SMesh, CGAL::vertex_point_t>::type Point_3;
    typedef typename boost::graph_traits<SMesh>::vertex_descriptor vertex_descriptor;

    VPMap vpmap = get(CGAL::vertex_point, tmesh);
    // get nb of points and cells
    vtkIdType nb_points = poly_data->GetNumberOfPoints();
    vtkIdType nb_cells = poly_data->GetNumberOfCells();
    //extract points
    std::vector<vertex_descriptor> vertex_map(nb_points);
    for (vtkIdType i = 0; i<nb_points; ++i)
    {
        double coords[3];
        poly_data->GetPoint(i, coords);

        vertex_descriptor v = CGAL::add_vertex(tmesh);
        put(vpmap, v, Point_3(coords[0], coords[1], coords[2]));
        vertex_map[i]=v;
    }
    for (vtkIdType i = 0; i<nb_cells; ++i)
    {
        if(poly_data->GetCellType(i) != 5
                && poly_data->GetCellType(i) != 7
                && poly_data->GetCellType(i) != 9) //only supported cells are triangles, quads and polygons
            continue;
        vtkCell* cell_ptr = poly_data->GetCell(i);

        vtkIdType nb_vertices = cell_ptr->GetNumberOfPoints();
        if (nb_vertices < 3)
            return false;
        std::vector<vertex_descriptor> vr(nb_vertices);
        for (vtkIdType k=0; k<nb_vertices; ++k)
            vr[k]=vertex_map[cell_ptr->GetPointId(k)];

        CGAL::Euler::add_face(vr, tmesh);
    }
    return true;
}

CGAL::Surface_mesh to vtkPolyData

首先，转换成vtkUnstructuredGrid格式，然后再转成vtkPolyData。

vtkSmartPointer<vtkUnstructuredGrid> polygon_mesh_to_vtkUnstructured(const SMesh &pmesh)
{
    typedef typename boost::graph_traits<SMesh>::vertex_descriptor   vertex_descriptor;
    typedef typename boost::graph_traits<SMesh>::face_descriptor     face_descriptor;
    typedef typename boost::graph_traits<SMesh>::halfedge_descriptor halfedge_descriptor;

    typedef typename boost::property_map<SMesh, CGAL::vertex_point_t>::const_type VPMap;
    typedef typename boost::property_map_value<SMesh, CGAL::vertex_point_t>::type Point_3;
    VPMap vpmap = get(CGAL::vertex_point, pmesh);

    vtkPoints* const vtk_points = vtkPoints::New();
    vtkCellArray* const vtk_cells = vtkCellArray::New();

    vtk_points->Allocate(CGAL::num_vertices(pmesh));
    vtk_cells->Allocate(CGAL::num_faces(pmesh));

    std::map<vertex_descriptor, vtkIdType> Vids;
    vtkIdType inum = 0;
    for(vertex_descriptor v : CGAL::vertices(pmesh))
    {
        const Point_3& p = get(vpmap, v);
        vtk_points->InsertNextPoint(CGAL::to_double(p.x()),
                                    CGAL::to_double(p.y()),
                                    CGAL::to_double(p.z()));
        Vids[v] = inum++;
    }
    for(face_descriptor f : CGAL::faces(pmesh))
    {
        vtkIdList* cell = vtkIdList::New();
        for(halfedge_descriptor h :
            CGAL::halfedges_around_face(CGAL::halfedge(f, pmesh), pmesh))
        {
            cell->InsertNextId(Vids[CGAL::target(h, pmesh)]);
        }
        vtk_cells->InsertNextCell(cell);
        cell->Delete();
    }

    vtkSmartPointer<vtkUnstructuredGrid> usg =
            vtkSmartPointer<vtkUnstructuredGrid>::New();

    usg->SetPoints(vtk_points);
    vtk_points->Delete();

    usg->SetCells(5,vtk_cells);
    vtk_cells->Delete();
    return usg;
}
// vtkUnstructuredGrid->vtkPolyData
vtkSmartPointer<vtkPolyData> ugToPolyData(vtkSmartPointer<vtkUnstructuredGrid> ug)
{
    vtkSmartPointer<vtkDataSetSurfaceFilter> surfaceFilter =
            vtkSmartPointer<vtkDataSetSurfaceFilter>::New();
    surfaceFilter->SetInputData(ug);
    surfaceFilter->Update();
    return surfaceFilter->GetOutput();
}

利用CGAL进行Mesh的布尔运算

//! op
//! enum bool_op {CRF_UNION(并), CRF_INTER（交）, CRF_MINUS（差）, CRF_MINUS_OP（差或）};
//！@return 布尔运算的结果Mesh

SMesh *apply_corefine_and_bool_op(SMesh *oneMesh, SMesh *twoMesh, bool_op op)
{
    if(! CGAL::is_triangle_mesh(*oneMesh)) {
        return nullptr;
    }
    if(! CGAL::is_triangle_mesh(*twoMesh)) {
        return nullptr;
    }
    SMesh* new_poly = new SMesh();
    QString str_op;
    SMesh P, Q;
    switch(op)
    {
    case CRF_UNION:
        P = *oneMesh, Q = *twoMesh;
        if (! PMP::corefine_and_compute_union(P, Q, *new_poly, params::throw_on_self_intersection(true)) )
        {
            delete new_poly;
            return nullptr;
        }
        str_op = "Union";
        break;
    case CRF_INTER:
        P = *oneMesh, Q = *twoMesh;
        if (! PMP::corefine_and_compute_intersection(P, Q, *new_poly, params::throw_on_self_intersection(true)) )
        {
            delete new_poly;
            return nullptr;
        }
        str_op = "Intersection";
        break;
    case CRF_MINUS_OP:
        std::swap(oneMesh, twoMesh);
        CGAL_FALLTHROUGH;
    case CRF_MINUS:
        P = *oneMesh, Q = *twoMesh;
        if (! PMP::corefine_and_compute_difference(P, Q, *new_poly, params::throw_on_self_intersection(true)) )
        {
            delete new_poly;
            return nullptr;
        }
        str_op = "Difference";
    }
    return new_poly;
}

Usage

    SMesh* meshOne = new SMesh();
    vtkPointSet_to_polygon_mesh(one,*meshOne);
    SMesh* meshTwo = new SMesh();
    vtkPointSet_to_polygon_mesh(two,*meshTwo);

    SMesh* meshThree=apply_corefine_and_bool_op(meshOne,meshTwo,op);
    vtkSmartPointer<vtkUnstructuredGrid>ug=polygon_mesh_to_vtkUnstructured(*meshThree);
    vtkSmartPointer<vtkPolyData>outPolyData= ugToPolyData(ug);
    delete meshOne;
    delete meshTwo;
    delete meshThree;

Result