CGAL users discussion list

["Laurent Rineau (GeometryFactory)" <>]

On Friday 12 March 2010 18:46:44 Ramin H wrote:
> Hi,
> 
> f.first returns a "cell handle" and f.second is the facet index in that
> cell. So f.first->vertex((f.second)&3) should be the vertex opposite to
> facet f, right?
> 
> Say p0 is the point opposite to the facet and p1, p2, p3 are on the facet,
> then:
> - We find "v1", the vector from p0 to p1
> - We find "n", the facet normal vector
> - To get a coherent direction for all normal vectors, we find the inner
> product of "n" and "v1" vectors. Only if the inner product is negative, we
> change the direction of the normal vector.
> 
> I have included the pseudocode below. Does this make sense?

> // ---- pseudocode ----
> // the point opposite to facet
> const Point_3& p0 = f.first->vertex((f.second)&3)->point();
> // points on the facet
> const Point_3& p1 = f.first->vertex((f.second+1)&3)->point();
> const Point_3& p2 = f.first->vertex((f.second+2)&3)->point();
> const Point_3& p3 = f.first->vertex((f.second+3)&3)->point();
> 
> Vector_3 v1 = p1 - p0;
> Vector_3 n = CGAL::normal(p1, p2, p3);
> 
> if ( n * v1 < 0 ) { n = -n; }
> 
> // ---> At this point we should have the correct direction for normal
> vector?

That seems right. The last test (n *v1 < 0) could be written:

  if(CGAL::angle(n, v1) == CGAL::ACUTE)

to be exact. But that would assume that n and v1 are constructed exactly, and 
that is false.

Actually, you could use the fact that cells are correctly oriented, in a CGAL 
3D triangulation. That way you can deduce the orientation of n without any 
numerical test. The parity of f.second is enough to know if the facet p1, p2, 
p3) is correctly oriented.

const Point_3& p0 = f.first->vertex((f.second)->point();
// points on the facet
const Point_3& p1 = f.first->vertex((f.second+1)&3)->point();
const Point_3& p2 = f.first->vertex((f.second+2)&3)->point();
const Point_3& p3 = f.first->vertex((f.second+3)&3)->point();

const Vector_3 n = ( f.second % 2 == 1) ? 
                   CGAL::normal(p1, p2, p3) :
                   CGAL::normal(p1, p3, p2);


You can check, in a positively oriented tetrahedron with vertices #0 #1 #2 
#3, 
then the positive facets are:
  #1 #2 #3
  #2 #0 #3
  #3 #0 #1
  #0 #2 #1

-- 
Laurent Rineau, PhD
R&D Engineer at GeometryFactory           http://www.geometryfactory.com/
Release Manager of the CGAL Project       http://www.cgal.org/