polymesht.cc

penMesh is a generic and efficient data structure for representing and manipulating polygonal meshes

//=============================================================================
//
//                               OpenMesh
//      Copyright (C) 2001-2005 by Computer Graphics Group, RWTH Aachen
//                           www.openmesh.org
//
//-----------------------------------------------------------------------------
//
//                                License
//
//   This library is free software; you can redistribute it and/or modify it
//   under the terms of the GNU Library General Public License as published
//   by the Free Software Foundation, version 2.
//
//   This library is distributed in the hope that it will be useful, but
//   WITHOUT ANY WARRANTY; without even the implied warranty of
//   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
//   Library General Public License for more details.
//
//   You should have received a copy of the GNU Library General Public
//   License along with this library; if not, write to the Free Software
//   Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
//
//-----------------------------------------------------------------------------
//
//   $Revision: 1.5 $
//   $Date: 2005-12-21 13:51:49 $
//
//=============================================================================


//=============================================================================
//
//  CLASS PolyMeshT - IMPLEMENTATION
//
//=============================================================================


#define OPENMESH_POLYMESH_C


//== INCLUDES =================================================================

#include <OpenMesh/Core/Mesh/PolyMeshT.hh>
#include <OpenMesh/Core/Math/LoopSchemeMaskT.hh>
#include <OpenMesh/Core/Utils/vector_cast.hh>
#include <OpenMesh/Core/System/omstream.hh>
#include <vector>


//== NAMESPACES ===============================================================


namespace OpenMesh {


//== IMPLEMENTATION ==========================================================


template <class Kernel>
const typename PolyMeshT<Kernel>::VertexHandle
PolyMeshT<Kernel>::InvalidVertexHandle;

template <class Kernel>
const typename PolyMeshT<Kernel>::HalfedgeHandle
PolyMeshT<Kernel>::InvalidHalfedgeHandle;

template <class Kernel>
const typename PolyMeshT<Kernel>::EdgeHandle
PolyMeshT<Kernel>::InvalidEdgeHandle;

template <class Kernel>
const typename PolyMeshT<Kernel>::FaceHandle
PolyMeshT<Kernel>::InvalidFaceHandle;



template <class Kernel_>
uint PolyMeshT<Kernel_>::find_feature_edges(Scalar _angle_tresh)
{
  assert(This::has_edge_status());//this function needs edge status property
  uint n_feature_edges = 0;
  for (EdgeIter e_it = edges_begin(); e_it != edges_end(); ++e_it)
  {
    if (fabs(calc_dihedral_angle(e_it)) > _angle_tresh)
    {//note: could be optimized by comparing cos(dih_angle) vs. cos(_angle_tresh)
      status(e_it).set_feature(true);
      n_feature_edges++;
    }
    else
    {
      status(e_it).set_feature(false);
    }
  }
  return n_feature_edges;
}

//-----------------------------------------------------------------------------


template <class Kernel_>
typename PolyMeshT<Kernel_>::HalfedgeHandle
PolyMeshT<Kernel_>::
find_halfedge( VertexHandle _start_vertex_handle,
               VertexHandle _end_vertex_handle ) const
{
  assert(_start_vertex_handle.is_valid() && _end_vertex_handle.is_valid());

  for (ConstVertexVertexIter vvIt=cvv_iter(_start_vertex_handle); vvIt; ++vvIt)
    if (vvIt.handle() == _end_vertex_handle)
      return vvIt.current_halfedge_handle();

  return InvalidHalfedgeHandle;
}


//-----------------------------------------------------------------------------


template <class Kernel>
void
PolyMeshT<Kernel>::
adjust_outgoing_halfedge(VertexHandle _vh)
{
  for (ConstVertexOHalfedgeIter vh_it=cvoh_iter(_vh); vh_it; ++vh_it)
  {
    if (is_boundary(vh_it.handle()))
    {
      set_halfedge_handle(_vh, vh_it.handle());
      break;
    }
  }
}


//-----------------------------------------------------------------------------


namespace add_face_ns
{
  template <class _Handle>
  struct NextCacheEntryT : public std::pair<_Handle, _Handle>
  {
    typedef std::pair<_Handle, _Handle> Base;

    NextCacheEntryT(_Handle _heh0, _Handle _heh1)
      : Base(_heh0, _heh1)
    {
      assert(_heh0.is_valid());
      assert(_heh1.is_valid());
    }
  };
};


template <class Kernel>
typename PolyMeshT<Kernel>::FaceHandle
PolyMeshT<Kernel>::
add_face(const std::vector<VertexHandle>& _vertex_handles)
{
  VertexHandle                   vh;
  unsigned int                   i, ii, n(_vertex_handles.size()), id;
  std::vector<HalfedgeHandle>    halfedge_handles(n);
  std::vector<bool>              is_new(n), needs_adjust(n, false);
  HalfedgeHandle                 inner_next, inner_prev,
                                 outer_next, outer_prev,
                                 boundary_next, boundary_prev,
                                 patch_start, patch_end;

  // cache for set_next_halfedge and vertex' set_halfedge
  typedef add_face_ns::NextCacheEntryT<HalfedgeHandle> NextCacheEntry;
  typedef std::vector<NextCacheEntry>                  NextCache;

  NextCache    next_cache;
  next_cache.reserve(3*n);


  // don't allow degenerated faces
  assert (n > 2);


  // test for topological errors
  for (i=0, ii=1; i<n; ++i, ++ii, ii%=n)
  {
    if ( !is_boundary(_vertex_handles[i]) )
    {
      omerr() << "PolyMeshT::add_face: complex vertex\n";
      return InvalidFaceHandle;
    }

    halfedge_handles[i] = find_halfedge(_vertex_handles[i],
                                        _vertex_handles[ii]);
    is_new[i] = !halfedge_handles[i].is_valid();

    if (!is_new[i] && !is_boundary(halfedge_handles[i]))
    {
      omerr() << "PolyMeshT::add_face: complex edge\n";
      return InvalidFaceHandle;
    }
  }


  // re-link patches if necessary
  for (i=0, ii=1; i<n; ++i, ++ii, ii%=n)
  {
    if (!is_new[i] && !is_new[ii])
    {
      inner_prev = halfedge_handles[i];
      inner_next = halfedge_handles[ii];

      if (next_halfedge_handle(inner_prev) != inner_next)
      {
        // here comes the ugly part... we have to relink a whole patch

        // search a free gap
        // free gap will be between boundary_prev and boundary_next
        outer_prev = opposite_halfedge_handle(inner_next);
        outer_next = opposite_halfedge_handle(inner_prev);
        boundary_prev = outer_prev;
        do
          boundary_prev =
            opposite_halfedge_handle(next_halfedge_handle(boundary_prev));
        while (!is_boundary(boundary_prev) || boundary_prev==inner_prev);
        boundary_next = next_halfedge_handle(boundary_prev);
        assert(is_boundary(boundary_prev));
        assert(is_boundary(boundary_next));


        // ok ?
        if (boundary_next == inner_next)
        {
          omerr() << "PolyMeshT::add_face: patch re-linking failed\n";
          return InvalidFaceHandle;
        }

        // other halfedges' handles
        patch_start = next_halfedge_handle(inner_prev);
        patch_end   = prev_halfedge_handle(inner_next);

        // relink
        next_cache.push_back(NextCacheEntry(boundary_prev, patch_start));
        next_cache.push_back(NextCacheEntry(patch_end, boundary_next));
        next_cache.push_back(NextCacheEntry(inner_prev, inner_next));
      }
    }
  }



  // create missing edges
  for (i=0, ii=1; i<n; ++i, ++ii, ii%=n)
    if (is_new[i])
      halfedge_handles[i] = new_edge(_vertex_handles[i], _vertex_handles[ii]);



  // create the face
  FaceHandle fh(This::new_face());
  set_halfedge_handle(fh, halfedge_handles[n-1]);



  // setup halfedges
  for (i=0, ii=1; i<n; ++i, ++ii, ii%=n)
  {
    vh         = _vertex_handles[ii];
    inner_prev = halfedge_handles[i];
    inner_next = halfedge_handles[ii];

    id = 0;
    if (is_new[i])  id |= 1;
    if (is_new[ii]) id |= 2;

    if (id)
    {
      outer_prev = opposite_halfedge_handle(inner_next);
      outer_next = opposite_halfedge_handle(inner_prev);

      // set outer links
      switch (id)
      {
        case 1: // prev is new, next is old
          boundary_prev = prev_halfedge_handle(inner_next);
          next_cache.push_back(NextCacheEntry(boundary_prev, outer_next));
          set_halfedge_handle(vh, outer_next);
          break;

        case 2: // next is new, prev is old
          boundary_next = next_halfedge_handle(inner_prev);
          next_cache.push_back(NextCacheEntry(outer_prev, boundary_next));
          set_halfedge_handle(vh, boundary_next);
          break;

        case 3: // both are new
          if (!halfedge_handle(vh).is_valid())
          {
            set_halfedge_handle(vh, outer_next);
            next_cache.push_back(NextCacheEntry(outer_prev, outer_next));
          }
          else
          {
            boundary_next = halfedge_handle(vh);
            boundary_prev = prev_halfedge_handle(boundary_next);
            next_cache.push_back(NextCacheEntry(boundary_prev, outer_next));
            next_cache.push_back(NextCacheEntry(outer_prev, boundary_next));
          }
          break;
      }

      // set inner link
      next_cache.push_back(NextCacheEntry(inner_prev, inner_next));
    }
    else needs_adjust[ii] = (halfedge_handle(vh) == inner_next);


    // set face handle
    set_face_handle(halfedge_handles[i], fh);
  }



  // process next halfedge cache
  typename NextCache::const_iterator
    ncIt(next_cache.begin()), ncEnd(next_cache.end());
  for (; ncIt != ncEnd; ++ncIt)
    set_next_halfedge_handle(ncIt->first, ncIt->second);



  // adjust vertices' halfedge handle
  for (i=0; i<n; ++i)
    if (needs_adjust[i])
      adjust_outgoing_halfedge(_vertex_handles[i]);



  return fh;
}


//-----------------------------------------------------------------------------


#if OM_OUT_OF_CLASS_TEMPLATE
template <typename Kernel>
template <typename OtherMesh>
PolyMeshT<Kernel>&
PolyMeshT<Kernel>::
assign(const OtherMesh& _rhs)
#  include "PolyMeshT_assign.hh"
#endif


//-----------------------------------------------------------------------------


template <class Kernel>
typename PolyMeshT<Kernel>::Normal
PolyMeshT<Kernel>::
calc_face_normal(FaceHandle _fh) const
{
  assert(halfedge_handle(_fh).is_valid());
  ConstFaceVertexIter fv_it(cfv_iter(_fh));

  const Point& p0(point(fv_it));  ++fv_it;
  const Point& p1(point(fv_it));  ++fv_it;
  const Point& p2(point(fv_it));

  return calc_face_normal(p0, p1, p2);
}


//-----------------------------------------------------------------------------


template <class Kernel>
typename PolyMeshT<Kernel>::Normal
PolyMeshT<Kernel>::
calc_face_normal(const Point& _p0,
		 const Point& _p1,
		 const Point& _p2) const
{
#if 1
  // The OpenSG <Vector>::operator -= () does not support the type Point
  // as rhs. Therefore use vector_cast at this point!!!
  // Note! OpenSG distinguishes between Normal and Point!!!
  Normal p1p0(_p0);  p1p0 -= vector_cast<Normal>(_p1);
  Normal p1p2(_p2);  p1p2 -= vector_cast<Normal>(_p1);

  Normal n    = cross(p1p2, p1p0);
  Scalar norm = n.length();

  // The expression ((n *= (1.0/norm)),n) is used because the OpenSG
  // vector class does not return self after component-wise
  // self-multiplication with a scalar!!!
  return (norm != Scalar(0)) ? ((n *= (Scalar(1)/norm)),n) : Normal(0,0,0);
#else
  Point p1p0 = _p0;  p1p0 -= _p1;
  Point p1p2 = _p2;  p1p2 -= _p1;

  Normal n = vector_cast<Normal>(cross(p1p2, p1p0));
  Scalar norm = n.length();

  return (norm != 0.0) ? n *= (1.0/norm) : Normal(0,0,0);
#endif
}


//-----------------------------------------------------------------------------


template <class Kernel>
void
PolyMeshT<Kernel>::
update_normals(VertexNormalMode _mode)
{
  if (This::has_face_normals())    update_face_normals();
  if (This::has_vertex_normals())  update_vertex_normals(_mode);
}


//-----------------------------------------------------------------------------


template <class Kernel>
void
PolyMeshT<Kernel>::
update_face_normals()
{
  FaceIter f_it(faces_begin()), f_end(faces_end());

  for (; f_it != f_end; ++f_it)
    set_normal(f_it.handle(), calc_face_normal(f_it.handle()));
}


//-----------------------------------------------------------------------------


template <class Kernel>
typename PolyMeshT<Kernel>::Normal
PolyMeshT<Kernel>::
calc_vertex_normal(VertexHandle _vh) const
{
  Normal n;
  calc_vertex_normal_fast(_vh,n);
  return n;
}


//-----------------------------------------------------------------------------


template <class Kernel_>
void PolyMeshT<Kernel_>::
calc_vertex_normal_fast(VertexHandle _vh, Normal& _n) const
{
  _n[0]=_n[1]=_n[2]=Scalar(0.0);
  for (ConstVertexFaceIter vf_it=cvf_iter(_vh); vf_it; ++vf_it)
    _n += normal(vf_it.handle());
  Scalar norm = _n.length();
  if (norm != 0.0f) _n *= (1.0f/norm);
}