CGAL带岛多边形三角化，并输出(*.ply)格式的模型
模型输出的关键是节点和索引
#include <CGAL/Triangulation_vertex_base_with_id_2.h>
#include <CGAL/Triangulation_face_base_with_info_2.h>
因此注意这两个泛型，对比不带信息的
#include <CGAL/Triangulation_vertex_base_2.h>
#include <CGAL/Triangulation_face_base_2.h>，这两个增加了部分信息作为拓展。
这样Vertex_handle就可以读取这部分拓展的信息。
心得：CGAL的泛型机制真的很强大，拓展性很好。
// AxModelDelaunay.cpp : 定义控制台应用程序的入口点。
// #include "stdafx.h"
#include "shapefil.h" #include "CGAL/exceptions.h"
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Constrained_Delaunay_triangulation_2.h>
#include <CGAL/Triangulation_vertex_base_with_id_2.h>
#include <CGAL/Triangulation_face_base_with_info_2.h>
#include <CGAL/Polygon_2.h>
#include <iostream>
struct FaceInfo2
{
FaceInfo2(){}
int nesting_level;
bool in_domain()
{
return nesting_level % == ;
}
}; typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Triangulation_vertex_base_with_id_2<K> Vb;
typedef CGAL::Triangulation_face_base_with_info_2<FaceInfo2, K> Fbb;
typedef CGAL::Constrained_triangulation_face_base_2<K, Fbb> Fb;
typedef CGAL::Triangulation_data_structure_2<Vb, Fb> TDS;
typedef CGAL::Exact_predicates_tag Itag;
typedef CGAL::Constrained_Delaunay_triangulation_2<K, TDS, Itag> CDT;
typedef CDT::Point Point;
typedef CGAL::Polygon_2<K> Polygon_2;
typedef CDT::Vertex_handle Vertex_handle; void mark_domains(CDT& ct, CDT::Face_handle start, int index, std::list<CDT::Edge>& border)
{
if (start->info().nesting_level != -){
return;
}
std::list<CDT::Face_handle> queue;
queue.push_back(start);
while (!queue.empty()){
CDT::Face_handle fh = queue.front();
queue.pop_front();
if (fh->info().nesting_level == -){
fh->info().nesting_level = index;
for (int i = ; i < ; i++){
CDT::Edge e(fh, i);
CDT::Face_handle n = fh->neighbor(i);
if (n->info().nesting_level == -){
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.
void mark_domains(CDT& cdt)
{
for (CDT::All_faces_iterator it = cdt.all_faces_begin(); it != cdt.all_faces_end(); ++it){
it->info().nesting_level = -;
}
std::list<CDT::Edge> border;
mark_domains(cdt, cdt.infinite_face(), , border);
while (!border.empty()){
CDT::Edge e = border.front();
border.pop_front();
CDT::Face_handle n = e.first->neighbor(e.second);
if (n->info().nesting_level == -){
mark_domains(cdt, n, e.first->info().nesting_level + , border);
}
}
} int _tmain(int argc, _TCHAR* argv[])
{
//读取shp
const char * pszShapeFile = "data\\walltest.shp";
SHPHandle hShp = SHPOpen(pszShapeFile, "r");
int nShapeType, nVertices;
int nEntities = ;
double* minB = new double[];
double* maxB = new double[];
SHPGetInfo(hShp, &nEntities, &nShapeType, minB, maxB);
printf("ShapeType:%d\n", nShapeType);
printf("Entities:%d\n", nEntities);
CDT cdt;
for (int i = ; i < nEntities; i++)
{
int iShape = i;
SHPObject *obj = SHPReadObject(hShp, iShape);
printf("--------------Feature:%d------------\n", iShape);
int parts = obj->nParts;
int* partStart = obj->panPartStart;
int verts = obj->nVertices;
printf("nParts:%d\n", parts);
printf("nVertices:%d\n", verts);
for (int k = ; k < parts; k++)
{
Polygon_2 polygon1;
int fromIdx = partStart[k];
int toIdx = fromIdx;
if (k<parts-)
{
toIdx = partStart[k + ];
}
else
{
toIdx = verts;
} for (size_t j = fromIdx; j < toIdx; j++)
{
double x = obj->padfX[j];
double y = obj->padfY[j];
//Point_2 pt(x, y);
printf("%f,%f;", x, y);
polygon1.push_back(Point(x, y)); }
cdt.insert_constraint(polygon1.vertices_begin(), polygon1.vertices_end(), true);
}
printf("\n");
} //construct two non-intersecting nested polygons
//Polygon_2 polygon1;
//polygon1.push_back(Point(0, 0));
//polygon1.push_back(Point(2, 0));
//polygon1.push_back(Point(2, 2));
//polygon1.push_back(Point(0, 2));
//Polygon_2 polygon2;
//polygon2.push_back(Point(0.5, 0.5));
//polygon2.push_back(Point(1.5, 0.5));
//polygon2.push_back(Point(1.5, 1.5));
//polygon2.push_back(Point(0.5, 1.5)); ////Insert the polygons into a constrained triangulation
//CDT cdt;
//cdt.insert_constraint(polygon1.vertices_begin(), polygon1.vertices_end(), true);
//cdt.insert_constraint(polygon2.vertices_begin(), polygon2.vertices_end(), true); //Mark facets that are inside the domain bounded by the polygon
mark_domains(cdt);
FILE *ply = fopen("data\\floorpeint.ply", "w"); int idx = ;
for (CDT::Vertex_iterator v = cdt.vertices_begin(); v != cdt.vertices_end(); ++v)
{
Vertex_handle vv = v->handle();
vv->id() = idx;
idx++;
} int count = ;
for (CDT::Finite_faces_iterator fit = cdt.finite_faces_begin();
fit != cdt.finite_faces_end(); ++fit)
{
if (fit->info().in_domain())
{
++count;
for (int i = ; i < ; i++)
{
Vertex_handle vert = fit->vertex(i);
int x=vert->id();
std::cout << "The Id is " << x << std::endl;
CDT::Edge ed(fit, i);
ed.second;
}
} }
if (ply)
{
fprintf(ply, "ply\nformat %s 1.0\n", "ascii");
fprintf(ply, "element vertex %d\n",idx );
fprintf(ply, "property float x\n");
fprintf(ply, "property float y\n");
fprintf(ply, "property float z\n");
fprintf(ply, "element face %d\n", count);
fprintf(ply, "property list uint8 int32 vertex_indices\n");
fprintf(ply, "end_header\n"); for (CDT::Vertex_iterator v = cdt.vertices_begin(); v != cdt.vertices_end(); ++v)
{
Vertex_handle vv = v->handle();
double x = vv->point().x();
double y = vv->point().y();
fprintf(ply, "%f %f %f\n", x, y, 0.0);
}
for (CDT::Finite_faces_iterator fit = cdt.finite_faces_begin(); fit != cdt.finite_faces_end(); ++fit)
{
if (fit->info().in_domain())
{ Vertex_handle vertId0 = fit->vertex();
Vertex_handle vertId1 = fit->vertex();
Vertex_handle vertId2 = fit->vertex();
int id0 = vertId0->id();
int id1 = vertId1->id();
int id2 = vertId2->id();
fprintf(ply, "%d %d %d %d\n", , id0, id1, id2);
} }
} std::cout << "There are " << count << " facets in the domain." << std::endl;
system("pause");
return ;
}
效果图