Hi all,

I want to be able to mesh 3d surfaces from an implicit function representing

several smooth closed surfaces.

From what I understand by reading the doc and a previous post

(

http://cgal-discuss.949826.n4.nabble.com/About-the-surface-mesh-generator-td952799.html)

I need to add one initial sample point for each connected component to the

initial triangulation.

I did a test with two spheres (inspired from CGAL tests

@

https://github.com/CGAL/cgal/blob/69fad29842b8a3a7316ea478b3a4da3a47d6bf6a/Surface_mesher/test/Surface_mesher/implicit_surface_mesher_test.cpp#L138),

but I can't get the final mesh with both the spheres in it.

<

http://cgal-discuss.949826.n4.nabble.com/file/t376135/Capture.png>

Is there something else I need to do after adding the initial points to the

triangulation or am I doing something wrong?

Thanks for your help,

What I did (see code below):

1 . create two implicit spheres

- left -> center = (-4,0,0), radius = 2

- right -> center = (4, 0,0), radius = 2

2. Insert one seed point per sphere to the triangulation

- left -> seed = (-2,0,0);

- right -> seed = (2,0,0);

3. Mesh the implicit surface representing both spheres

4. Output mesh to .off file

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>

#include <CGAL/Surface_mesh_default_triangulation_3.h>

#include <CGAL/Complex_2_in_triangulation_3.h>

#include <CGAL/make_surface_mesh.h>

#include <CGAL/Implicit_surface_3.h>

#include <CGAL/Polyhedron_3.h>

#include <CGAL/IO/facets_in_complex_2_to_triangle_mesh.h>

#include <string>

#include <iostream>

#include <fstream>

template <typename K>

class Sphere {

public:

typedef typename K::FT FT;

typedef typename K::Point_3 Point_3;

typedef typename K::Sphere_3 Sphere_3;

Sphere(Sphere_3 sphere)

: sphere(sphere)

{

}

FT operator()(const Point_3& p) const

{

FT x = p.x();

FT y = p.y();

FT z = p.z();

x-=sphere.center().x();

y-=sphere.center().y();

z-=sphere.center().z();

return x*x+y*y+z*z-sphere.squared_radius();

}

private:

const Sphere_3 sphere;

};

template <typename K>

class TwoSpheres

{

public:

typedef typename K::FT FT;

typedef typename K::Point_3 Point_3;

TwoSpheres(Sphere<K> sphere1, Sphere<K> sphere2)

: sphere1(sphere1), sphere2(sphere2)

{

}

FT operator()(const Point_3& p) const

{

return sphere1(p)*sphere2(p);

}

private:

const Sphere<K> sphere1;

const Sphere<K> sphere2;

};

// Aliases

using Kernel =

CGAL::Exact_predicates_inexact_constructions_kernel;

using Sphere_3 = Kernel::Sphere_3;

using Point_3 = Kernel::Point_3;

using FT = Kernel::FT;

using Tr = CGAL::Surface_mesh_default_triangulation_3;

using C2t3 = CGAL::Complex_2_in_triangulation_3

;

using Polyhedron = CGAL::Polyhedron_3<Kernel>;

using ImplicitSurface_3 = CGAL::Implicit_surface_3<Kernel,

TwoSpheres<Kernel>>;

int main() {

// Construct the spheres

// - left -> center = (-4,0,0), radius = 2

// - right -> center = (4, 0,0), radius = 2

const Sphere<Kernel> sphere_left (Sphere_3(Point_3(-4., 0., 0.), 4.));

const Sphere<Kernel> sphere_right(Sphere_3(Point_3(4., 0., 0.), 4.));

// Assemble the two spheres

const TwoSpheres<Kernel> two_spheres(sphere_left, sphere_right);

// Create the seeds for the 3D-Delaunay triangulation

// - one seed per connected component

const Point_3 point_left (-2,0,0);

const Point_3 point_right(2,0,0);

if(two_spheres(point_left) != 0)

std::cerr << "ERROR: Left seed is not on sphere!" << std::endl;

if(two_spheres(point_right) != 0)

std::cerr << "ERROR: Right seed is not on sphere!" << std::endl;

// Set the 3D-Delaunay triangulation

Tr tr;

tr.insert(point_left);

tr.insert(point_right);

// Set the 2D-complex in 3D-Delaunay triangulation

C2t3 c2t3 (tr);

// Define the implicit surface

const ImplicitSurface_3 surface(two_spheres,

Sphere_3(CGAL::ORIGIN,

100));

// Define meshing criteria

const CGAL::Surface_mesh_default_criteria_3

criteria(30., // angular bound

0.1, // radius bound

0.1); // distance bound

// Mesh surface

CGAL::make_surface_mesh(c2t3, surface, criteria,

CGAL::Non_manifold_tag());

std::cout << "Final number of points: " << tr.number_of_vertices() <<

std::endl;

// Converts to polyhedron

Polyhedron output_mesh;

CGAL::facets_in_complex_2_to_triangle_mesh(c2t3, output_mesh);

// Save to .off file

const std::string output_filename = "test.off";

std::cout << "Write file " << output_filename << std::endl;

std::ofstream out(output_filename.c_str());

out << output_mesh;

}

--

Sent from:

http://cgal-discuss.949826.n4.nabble.com/--

You are currently subscribed to cgal-discuss.

To unsubscribe or access the archives, go to

https://sympa.inria.fr/sympa/info/cgal-discuss