sqrt3t.hh

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.3.2.1 $
//   $Date: 2007-02-26 12:51:46 $
//                                                                            
//=============================================================================

/** \file Sqrt3T.hh
    
 */

//=============================================================================
//
//  CLASS Sqrt3T
//
//=============================================================================

#ifndef OPENMESH_SUBDIVIDER_UNIFORM_SQRT3T_HH
#define OPENMESH_SUBDIVIDER_UNIFORM_SQRT3T_HH


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

#include <OpenMesh/Core/System/config.hh>
#include <OpenMesh/Tools/Subdivider/Uniform/SubdividerT.hh>
#if defined(_DEBUG) || defined(DEBUG)
// Makes life lot easier, when playing/messing around with low-level topology
// changing methods of OpenMesh
#  include <OpenMesh/Tools/Utils/MeshCheckerT.hh>
#  define ASSERT_CONSISTENCY( T, m ) \
     assert(OpenMesh::Utils::MeshCheckerT<T>(m).check())
#else
#  define ASSERT_CONSISTENCY( T, m )
#endif
// -------------------- STL
#include <vector>
#if defined(OM_CC_MIPS)
#  include <math.h>
#else
#  include <cmath>
#endif


//== NAMESPACE ================================================================

namespace OpenMesh   { // BEGIN_NS_OPENMESH
namespace Subdivider { // BEGIN_NS_DECIMATER
namespace Uniform    { // BEGIN_NS_DECIMATER


//== CLASS DEFINITION =========================================================


/** %Uniform Sqrt3 subdivision algorithm
 *
 *  Implementation as described in
 *
 *  L. Kobbelt, <a href="http://www-i8.informatik.rwth-aachen.de/publications/downloads/sqrt3.pdf">"Sqrt(3) subdivision"</a>, Proceedings of SIGGRAPH 2000.
 */
template <typename MeshType, typename RealType = float>
class Sqrt3T : public SubdividerT< MeshType, RealType >
{
public:

  typedef RealType                                real_t;
  typedef MeshType                                mesh_t;
  typedef SubdividerT< mesh_t, real_t >           parent_t;

  typedef std::pair< real_t, real_t >             weight_t;
  typedef std::vector< std::pair<real_t,real_t> > weights_t;

public:


  Sqrt3T(void) : parent_t(), _1over3( 1.0/3.0 ), _1over27( 1.0/27.0 )
  { init_weights(); }

  virtual ~Sqrt3T() {}


public:


  const char *name() const { return "Uniform Sqrt3"; }

  
  /// Pre-compute weights
  void init_weights(size_t _max_valence=50)
  {
    weights_.resize(_max_valence);
    std::generate(weights_.begin(), weights_.end(), compute_weight());
  }


protected:


  bool prepare( MeshType& _m )
  {
    _m.request_edge_status();
    _m.add_property( vp_pos_ );
    _m.add_property( ep_nv_ );
    _m.add_property( mp_gen_ );
    _m.property( mp_gen_ ) = 0;

    return _m.has_edge_status() && vp_pos_.is_valid() 
      &&   ep_nv_.is_valid() && mp_gen_.is_valid();
  }


  bool cleanup( MeshType& _m )
  {
    _m.release_edge_status();
    _m.remove_property( vp_pos_ );
    _m.remove_property( ep_nv_ );
    _m.add_property( mp_gen_ );
    return true;
  }


  bool subdivide( MeshType& _m, size_t _n )
  {
    typename MeshType::VertexIter       vit;
    typename MeshType::VertexVertexIter vvit;
    typename MeshType::EdgeIter         eit;
    typename MeshType::FaceIter         fit;
    typename MeshType::FaceVertexIter   fvit;
    typename MeshType::VertexHandle     vh;
    typename MeshType::HalfedgeHandle   heh;
    typename MeshType::Point            pos(0,0,0), zero(0,0,0);
    size_t                            &gen = _m.property( mp_gen_ );

    for (size_t l=0; l<_n; ++l)
    {
      // tag existing edges
      for (eit=_m.edges_begin(); eit != _m.edges_end();++eit)
      {
        _m.status( eit ).set_tagged( true );
        if ( (gen%2) && _m.is_boundary(eit) )
          compute_new_boundary_points( _m, eit ); // *) creates new vertices
      }

      // do relaxation of old vertices, but store new pos in property vp_pos_

      for (vit=_m.vertices_begin(); vit!=_m.vertices_end(); ++vit)
      {
        if ( _m.is_boundary(vit) )
        {
          if ( gen%2 )
          {
            heh  = _m.halfedge_handle(vit);
            if (heh.is_valid()) // skip isolated newly inserted vertices *)
            {
              typename MeshType::HalfedgeHandle 
                prev_heh = _m.prev_halfedge_handle(heh);

              assert( _m.is_boundary(heh     ) );
              assert( _m.is_boundary(prev_heh) );
            
              pos  = _m.point(_m.to_vertex_handle(heh));
              pos += _m.point(_m.from_vertex_handle(prev_heh));
              pos *= real_t(4.0);

              pos += real_t(19.0) * _m.point( vit );
              pos *= _1over27;

              _m.property( vp_pos_, vit ) = pos;
            }
          }
          else
            _m.property( vp_pos_, vit ) = _m.point( vit );
        }
        else
        {
          size_t valence=0;

          pos = zero;
          for ( vvit = _m.vv_iter(vit); vvit; ++vvit)
          {
            pos += _m.point( vvit );
            ++valence;
          }
          pos *= weights_[ valence ].second;
          pos += weights_[ valence ].first * _m.point(vit);
          _m.property( vp_pos_, vit ) =  pos;
        }
      }   

      // insert new vertices, but store pos in vp_pos_
      typename MeshType::FaceIter fend = _m.faces_end();
      for (fit = _m.faces_begin();fit != fend; ++fit)
      {
        if ( (gen%2) && _m.is_boundary(fit))
        {
          boundary_split( _m, fit );
        }
        else
        {
          fvit = _m.fv_iter( fit );        
          pos  = _m.point(  fvit);
          pos += _m.point(++fvit);
          pos += _m.point(++fvit);
          pos *= _1over3;
          vh   = _m.add_vertex( zero );
          _m.property( vp_pos_, vh ) = pos;
          _m.split( fit, vh );
        }
      }

      // commit new positions (now iterating over all vertices)
      for (vit=_m.vertices_begin();vit != _m.vertices_end(); ++vit)
        _m.set_point(vit, _m.property( vp_pos_, vit ) );
      
      // flip old edges
      for (eit=_m.edges_begin(); eit != _m.edges_end(); ++eit)
        if ( _m.status( eit ).tagged() && !_m.is_boundary( eit ) )
          _m.flip(eit);

      // Now we have an consistent mesh!
      ASSERT_CONSISTENCY( MeshType, _m );

      // increase generation by one
      ++gen;
    }
    return true;
  }

private: