optimize triangulate and fix normal

.gitignore

@@ -18,3 +18,5 @@ polygon_triangulation2.obj
 
 gmon.out
 
+include/SFCGAL/version.h
+

CMakeLists.txt

@@ -6,6 +6,7 @@ set( CMAKE_DEBUG_POSTFIX "d" )
 #-- include finders and co
 set( CMAKE_MODULE_PATH "${CMAKE_SOURCE_DIR}/cmake/Modules;${CMAKE_MODULE_PATH}" )
 
+include( sfcgal-config.cmake )
 include( PrecompiledHeader )
 option( Use_precompiled_headers "Use precompiled headers" OFF )
 
@@ -105,8 +106,9 @@ elseif(CMAKE_COMPILER_IS_GNUCC OR CMAKE_COMPILER_IS_GNUCXX)
 	set(CMAKE_EXE_LINKER_FLAGS_DEBUG "-pg")
 endif()
 
-#-- configure the library
+#-- generate library headers
 configure_file( ${CMAKE_SOURCE_DIR}/include/SFCGAL/config.h.cmake ${CMAKE_SOURCE_DIR}/include/SFCGAL/config.h )
+configure_file( ${CMAKE_SOURCE_DIR}/include/SFCGAL/version.h.cmake ${CMAKE_SOURCE_DIR}/include/SFCGAL/version.h )
 
 #-- build the library
 #note : not available on windows without export/import

include/SFCGAL/Polygon.h

@@ -36,6 +36,10 @@ namespace SFCGAL {
 		 * Constructor with an exterior ring
 		 */
 		Polygon( const LineString & exteriorRing ) ;
+		/**
+		 * Constructor with an exterior ring (takes ownership)
+		 */
+		Polygon( LineString * exteriorRing ) ;
 		/**
 		 * Constructor with a Triangle
 		 */
@@ -97,7 +101,18 @@ namespace SFCGAL {
 		inline LineString &          exteriorRing() {
 			return _rings.front();
 		}
-
+		/**
+		 * Sets the exterior ring
+		 */
+		inline void  setExteriorRing( const LineString& ring ){
+			_rings.front() = ring ;
+		}
+		/**
+		 * Sets the exterior ring (takes ownership)
+		 */
+		inline void  setExteriorRing( LineString* ring ){
+			_rings.replace( 0, ring );
+		}
 
 		/**
 		 * Test if the polygon has interior rings

include/SFCGAL/algorithm/normal.h

@@ -23,33 +23,39 @@ namespace algorithm {
 
 	/**
 	 * Returns the 3D normal to a ring (supposed to be planar and closed).
-	 */
+         * @warning exact allows to avoid double rounding at the end of the computation
+	*/
 	template < typename Kernel >
-	CGAL::Vector_3< Kernel > normal3D( const LineString & ls )
+	CGAL::Vector_3< Kernel > normal3D( const LineString & ls, bool exact = true )
 	{
 		// Newell's formula
 		typename Kernel::FT nx, ny, nz;
 		nx = ny = nz = 0.0;
-		for ( size_t i = 0; i < ls.numPoints() - 1; ++i )
+		for ( size_t i = 0; i < ls.numPoints(); ++i )
 		{
 			const Point& pi = ls.pointN(i);
 			const Point& pj = ls.pointN( (i+1) % ls.numPoints() );
-			typename Kernel::FT zi = pi.is3D() ? pi.z() : 0.0;
-			typename Kernel::FT zj = pj.is3D() ? pj.z() : 0.0;
+			typename Kernel::FT zi = pi.z() ;
+			typename Kernel::FT zj = pj.z() ;
 			nx += ( pi.y() - pj.y() ) * ( zi + zj );
 			ny += ( zi - zj ) * ( pi.x() + pj.x() );
 			nz += ( pi.x() - pj.x() ) * ( pi.y() + pj.y() );
 		}
-		return CGAL::Vector_3<Kernel>( nx, ny, nz );
-	}
+		if ( exact ){
+			return CGAL::Vector_3<Kernel>( nx, ny, nz );
+		}else{
+			return CGAL::Vector_3<Kernel>( CGAL::to_double(nx), CGAL::to_double(ny), CGAL::to_double(nz) );
+		}	
+	} 
 
 	/**
 	 * Returns the 3D normal to a polygon (supposed to be planar).
+         * @warning exact allows to avoid double rounding at the end of the computation
 	 */
 	template < typename Kernel >
-	CGAL::Vector_3< Kernel > normal3D( const Polygon & polygon )
+	CGAL::Vector_3< Kernel > normal3D( const Polygon & polygon, bool exact = true )
 	{
-		return normal3D< Kernel >( polygon.exteriorRing() );
+		return normal3D< Kernel >( polygon.exteriorRing(), exact );
 	}
 
 }//algorithm

include/SFCGAL/algorithm/plane.h

@@ -40,10 +40,11 @@ namespace algorithm {
 			}else if ( n == 2 && ! CGAL::collinear( a, b, p ) ) {
 				c = p ;
 				n++ ;
-				break;
+				return true ;
 			}
 		}
-		return n == 3;
+		BOOST_ASSERT( n < 3 );
+		return false;
 	}
 
 	/**
@@ -79,18 +80,17 @@ namespace algorithm {
 
 	/**
 	 * Returns the oriented 3D plane of a polygon (supposed to be planar).
+	 * @warning result is rounded to double if exact is false (avoid huge expression tree)
 	 */
 	template < typename Kernel >
-	CGAL::Plane_3< Kernel > plane3D( const Polygon & polygon )
+	CGAL::Plane_3< Kernel > plane3D( const Polygon & polygon, bool exact = true )
 	{
-		CGAL::Vector_3< Kernel > normal = normal3D< Kernel >( polygon );
-		return CGAL::Plane_3< Kernel >( polygon.exteriorRing().pointN(0).toPoint_3(),
-						normal );
+		CGAL::Vector_3< Kernel > normal = normal3D< Kernel >( polygon, exact );
+		return CGAL::Plane_3< Kernel >( polygon.exteriorRing().pointN(0).toPoint_3(), normal );
 	}
 
 
 
-
 }//algorithm
 }//SFCGAL
 

include/SFCGAL/generator/disc.h

@@ -0,0 +1,26 @@
+#ifndef _SFCGAL_GENERATOR_DISC_H_
+#define _SFCGAL_GENERATOR_DISC_H_
+
+#include <memory>
+
+namespace SFCGAL {
+	class Point ;
+	class Polygon ;
+}
+
+namespace SFCGAL {
+namespace generator {
+
+	/**
+	 * Generate a discrete circle
+	 */
+	std::auto_ptr< Polygon > disc(
+		const Point& center,
+		const double & radius,
+		const unsigned int & nQuadrantSegments = 8U
+	) ;
+
+} // namespace generator
+} // namespace SFCGAL
+
+#endif

include/SFCGAL/generator/hoch.h

@@ -0,0 +1,21 @@
+#ifndef _SFCGAL_GENERATOR_HOCH_H_
+#define _SFCGAL_GENERATOR_HOCH_H_
+
+#include <memory>
+
+namespace SFCGAL {
+	class Polygon ;
+}
+
+namespace SFCGAL {
+namespace generator {
+
+	/**
+	 * generate hoch snowflake
+	 */
+	std::auto_ptr< Polygon > hoch( const unsigned int & order ) ;
+
+} // namespace generator
+} // namespace SFCGAL
+
+#endif

include/SFCGAL/version.h.cmake

@@ -0,0 +1,11 @@
+#ifndef _SFCGAL_VERSION_H_
+#define _SFCGAL_VERSION_H_
+
+#define SFCAL_VERSION_MAJOR @SFCGAL_VERSION_MAJOR@
+#define SFCAL_VERSION_MINOR @SFCGAL_VERSION_MINOR@
+#define SFCAL_VERSION_PATCH @SFCGAL_VERSION_PATCH@
+
+#define SFCAL_VERSION "@SFCGAL_VERSION@"
+
+#endif
+

sfcgal-config.cmake

@@ -0,0 +1,6 @@
+set( SFCGAL_VERSION_MAJOR 0 )
+set( SFCGAL_VERSION_MINOR 1 )
+set( SFCGAL_VERSION_PATCH 0 )
+
+set( SFCGAL_VERSION "${SFCGAL_VERSION_MAJOR}.${SFCGAL_VERSION_MINOR}.${SFCGAL_VERSION_PATCH}" )
+

src/SFCGAL/Polygon.cpp

@@ -39,6 +39,15 @@ Polygon::Polygon( const LineString & exteriorRing ):
 	_rings.push_back( exteriorRing.clone() );
 }
 
+///
+///
+///
+Polygon::Polygon( LineString * exteriorRing ):
+	Surface()
+{
+	_rings.push_back( exteriorRing );
+}
+
 ///
 ///
 ///

src/SFCGAL/algorithm/triangulate_polygon.cpp

@@ -210,15 +210,10 @@ void triangulate( const MultiPoint & geometry, TriangulatedSurface & triangulate
 	for ( size_t j = 0; j < geometry.numGeometries(); j++ ) {
 		const Point & point = geometry.pointN( j );
 
-		CGAL::Point_2< Kernel > p2d(
-			point.x(),
-			point.y()
-		);
-
 		/*
 		 * insert into triangulation
 		 */
-		CDT::Vertex_handle vh = cdt.insert( p2d );
+		CDT::Vertex_handle vh = cdt.insert( point.toPoint_2() );
 		vh->info().original = point ;
 	}
 
@@ -267,9 +262,9 @@ void triangulate( const Polygon & polygon, TriangulatedSurface & triangulatedSur
 
 //	std::cout << "---------------------------------------------------------" << std::endl ;
 //	std::cout << "triangulate polygon : " << polygon.asText() << std::endl;
-	CGAL::Plane_3< Kernel > polygonPlane = plane3D< Kernel >( polygon ) ;
+	CGAL::Plane_3< Kernel > polygonPlane = plane3D< Kernel >( polygon, false ) ;
 	if ( polygonPlane.is_degenerate() ){
-		BOOST_THROW_EXCEPTION( Exception( 
+		BOOST_THROW_EXCEPTION( Exception(
 			( boost::format( "can't find plane for polygon %s" ) % polygon.asText() ).str()
 		) );
 	}
@@ -283,7 +278,6 @@ void triangulate( const Polygon & polygon, TriangulatedSurface & triangulatedSur
 		CDT::Vertex_handle v_prev ;
 		for ( size_t j = 0; j < ring.numPoints(); j++ ) {
 			const Point & point = ring.pointN( j );
-			SFCGAL_DEBUG( boost::format( "insert point %s" ) % point.asText() );
 			CGAL::Point_3< Kernel > p3d = point.toPoint_3();
 
 			/*
@@ -319,7 +313,6 @@ void triangulate( const Polygon & polygon, TriangulatedSurface & triangulatedSur
 		if ( ! it->info().in_domain() ){
 			continue ;
 		}
-//		assert( it->is_valid() );
 
 		const Point & a = it->vertex(0)->info().original ;
 		const Point & b = it->vertex(1)->info().original ;
@@ -362,15 +355,10 @@ void triangulate2D( const Polygon & polygon, TriangulatedSurface & triangulatedS
 		for ( size_t j = 0; j < ring.numPoints(); j++ ) {
 			const Point & point = ring.pointN( j );
 
-			CGAL::Point_2< Kernel > p2d(
-				point.x(),
-				point.y()
-			);
-
 			/*
 			 * insert into triangulation
 			 */
-			CDT::Vertex_handle vh = cdt.insert( p2d );
+			CDT::Vertex_handle vh = cdt.insert( point.toPoint_2() );
 			vh->info().original = point ;
 
 			// filter first point

src/SFCGAL/generator/disc.cpp

@@ -0,0 +1,32 @@
+#include <SFCGAL/generator/disc.h>
+#include <SFCGAL/all.h>
+
+namespace SFCGAL {
+namespace generator {
+
+///
+///
+///
+std::auto_ptr< Polygon > disc(
+	const Point& center,
+	const double & radius,
+	const unsigned int & nQuadrantSegments
+)
+{
+	BOOST_ASSERT( nQuadrantSegments > 1 );
+
+	std::auto_ptr< LineString > exteriorRing( new LineString() ) ;
+
+	double dTheta = M_PI_4 / nQuadrantSegments ;
+	for ( size_t i = 0; i < nQuadrantSegments * 4; i++ ){
+		Kernel::Vector_2 p = center.toVector_2() + Kernel::Vector_2( cos(i*dTheta), sin(i*dTheta) ) ;
+		exteriorRing->addPoint( new Point( p.x(), p.y() ) ) ;
+	}
+	exteriorRing->addPoint( exteriorRing->startPoint() ) ;
+
+	return std::auto_ptr< Polygon >( new Polygon( exteriorRing.release() ) );
+}
+
+} // namespace generator
+} // namespace SFCGAL
+

src/SFCGAL/generator/hoch.cpp

@@ -0,0 +1,56 @@
+#include <SFCGAL/generator/hoch.h>
+
+#include <SFCGAL/all.h>
+
+namespace SFCGAL {
+namespace generator {
+
+std::vector< Kernel::Vector_2 > _hoch( const std::vector< Kernel::Vector_2 > & points ){
+	std::vector< Kernel::Vector_2 > result ;
+	//TODO result.reserve
+	size_t numPoints = points.size() ;
+	for ( size_t i = 0; i < numPoints; i++ ){
+		const Kernel::Vector_2& a = points[ i ] ;
+		const Kernel::Vector_2& b = points[ (i + 1) % numPoints ] ;
+
+		Kernel::Vector_2 ab = b - a ;
+		Kernel::Vector_2 normal( -ab.y(), ab.x() );
+
+		result.push_back( a );
+		result.push_back( a + ab / 3 );
+		result.push_back( a + ab / 2 + normal * sqrt(3.0) / 6.0 );
+		result.push_back( a + (ab * 2)/3 );
+	}
+	return result ;
+}
+
+///
+///
+///
+std::auto_ptr< Polygon > hoch( const unsigned int & order )
+{
+	std::vector< Kernel::Vector_2 > points ;
+	points.push_back( Kernel::Vector_2( 1.0, sqrt(3.0) ) );
+	points.push_back( Kernel::Vector_2( 2.0, 0.0 ) );
+	points.push_back( Kernel::Vector_2( 0.0, 0.0 ) );
+
+
+	for ( unsigned int k = 0 ; k < order; k++ ){
+		points = _hoch( points );
+	}
+
+	std::auto_ptr< Polygon > result( new Polygon() ) ;
+	std::auto_ptr< LineString > ring( new LineString() ) ;
+	for ( std::vector< Kernel::Vector_2 >::const_iterator it = points.begin(); it != points.end(); ++it ){
+		ring->addPoint( new Point( it->x(), it->y() ) ) ;
+	}
+	ring->addPoint( new Point( points.front().x(), points.front().y() ) ) ;
+
+	result->setExteriorRing( ring.release() );
+
+	return result ;
+}
+
+} // namespace generator
+} // namespace SFCGAL
+

test/CMakeLists.txt

@@ -12,8 +12,8 @@ add_subdirectory( unit )
 #-- build regression tests
 add_subdirectory( regress )
 
-option( SFCGAL_COMPILE_BENCH "Compile benchmarks" OFF )
-if( SFCGAL_COMPILE_BENCH )
+option( SFCGAL_BUILD_BENCH "Build benchmarks" ON )
+if( ${SFCGAL_BUILD_BENCH} )
 	add_subdirectory( bench )
 endif()