Computing normal of all facets in CGAL::Polyhedron_3

In this post I will show how to compute normal vector of a facet in CGAL::Polyhedron_3.

ComputeFacetNormal() Functor

ComputeFacetNormal() is a thread-safe functor for computing the normal vector of a given facet in a CGAL::Polyhedron_3. The facet must be a planar polygon with arbitrary number of sides.

#ifndef _SMESHLIB_OPERATIONS_COMPUTEFACETNORMAL_H_
#define _SMESHLIB_OPERATIONS_COMPUTEFACETNORMAL_H_

namespace SMeshLib   {
namespace Operations {
;

// The ComputeFacetNormal is a functor (delegate) to compute the normal of a facet in CGAL::Polyhdeon_3.
// ComputeFacetNormal is thread-safe.
// TPolyhedron is a type of CGAL::Polyhdeon_3.
template<class TPolyhedron>
struct ComputeFacetNormal
{
public:
	
	// Redefine types from TPolyhedron for convenience.
	typedef typename TPolyhedron::Facet                 Facet;
	typedef typename TPolyhedron::Halfedge_const_handle HalfedgeConstHandle;
	typedef typename TPolyhedron::Traits::Point_3       Point3;
	typedef typename TPolyhedron::Traits::Vector_3      Vector3;
	
	// Return type of operator() required by QtConcurrent.
	typedef Vector3 result_type;
	
public:
	
	// Compute normal of the given facet.
	// Facet can be triangle, quadrilateral or a polygon as long as its planar.
	// Use first three vertices to compute the normal.
	inline Vector3 operator() (const Facet& f) const
	{
		HalfedgeConstHandle h  = f.halfedge();
		Point3              p1 = h->vertex()->point();
		Point3              p2 = h->next()->vertex()->point();
		Point3              p3 = h->next()->next()->vertex()->point();
		Vector3             n  = CGAL::cross_product(p2-p1, p3-p1);
		return n / std::sqrt(n*n);
	}
};

};	// End namespace Operations.
};	// End namespace SMeshLib.

#endif // _SMESHLIB_OPERATIONS_COMPUTEFACETNORMAL_H_

Using ComputeFacetNormal() Functor

Normal vector of a facet f can be computed as Vector3 normal = ComputeFacetNormal(f);.

For most purposes, it is better to compute area of all facets once and cache them for later use. It is best to store the results in an associative container which associates the facet handle with the area. In the following example, I use PropertyMap which is a wrapper for std::set.

#include "ImportOBJ.h"
#include "ComputeFacetNormal.h"
#include "PropertyMap.h"
#include "CGAL/Simple_cartesian.h"
#include "CGAL/Polyhedron_items_3.h"
#include "CGAL/HalfedgeDS_list.h"
#include "CGAL/Polyhedron_3.h"

typedef CGAL::Simple_cartesian<double> Kernel;
typedef CGAL::Polyhedron_3<Kernel, 
	                   CGAL::Polyhedron_items_3, 
			   CGAL::HalfedgeDS_list> CgalPolyhedron;

// Compute the normal of all facets in a CGAL::Polyhedron_3 and stores it in 
// facetNormals.
// TPolyhedron is a type of CGAL::Polyhedron_3.
// TFacetNormals is a property map which associates a normal vector to facet in 
// the CGAL::Polyhdeon_3. 
template<class TPolyhedron, class TFacetNormals>
void computeFacetNormals(const TPolyhedron& polyhedron, TFacetNormals* facetNormals)
{
	if(facetNormals == 0)
	{
		return;
	}
	
	typename TPolyhedron::Facet_const_iterator _begin = polyhedron.facets_begin();
	typename TPolyhedron::Facet_const_iterator _end   = polyhedron.facets_end();

	SMeshLib::Operations::ComputeFacetNormal<TPolyhedron> _computeFacetNormal;
	CGAL_For_all(_begin, _end)
	{
		facetNormals->setValue(&*_begin, _computeFacetNormal(*_begin));
	}
}

void testComputeFacetNormal()
{
	CgalPolyhedron _poly;
	SMeshLib::IO::importOBJ("Venus.obj", &_poly);
	
	typedef SMeshLib::PropertyMap<const CgalPolyhedron::Facet*, CgalPolyhedron::Traits::Vector_3> FacetNormalPM;
	FacetNormalPM _facetNormals;
	computeFacetNormals(_poly, &_facetNormals);
}

Downloads

ImportOBJ.h
PropertyMap.h
ComputeFacetNormal.h
TestComputeFacetNormal.cpp
Venus.obj
ComputeFacetNormal.zip