sqrt3t.hh

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

  /// Helper functor to compute weights for sqrt(3)-subdivision
  /// \internal
  struct compute_weight 
  {
    compute_weight() : valence(-1) { }    
    weight_t operator() (void) 
    { 
#if !defined(OM_CC_MIPS)
      using std::cos;
#endif
      if (++valence)
      {
        real_t alpha = (4.0-2.0*cos(2.0*M_PI / (double)valence))/9.0;
        return weight_t( real_t(1)-alpha, alpha/real_t(valence) );
      }
      return weight_t(0.0, 0.0);
    }    
    int valence;
  };

private:

  // Pre-compute location of new boundary points for odd generations
  // and store them in the edge property ep_nv_;
  void compute_new_boundary_points( MeshType& _m, 
                                    const typename MeshType::EdgeHandle& _eh)
  {
    assert( _m.is_boundary(_eh) );

    typename MeshType::HalfedgeHandle heh;
    typename MeshType::VertexHandle   vh1, vh2, vh3, vh4, vhl, vhr;
    typename MeshType::Point          zero(0,0,0), P1, P2, P3, P4;

    /*
    //       *---------*---------*
    //      / \       / \       / \
    //     /   \     /   \     /   \
    //    /     \   /     \   /     \
    //   /       \ /       \ /       \
    //  *---------*--#---#--*---------*
    //                
    //  ^         ^  ^   ^  ^         ^
    //  P1        P2 pl  pr P3        P4
    */
    // get halfedge pointing from P3 to P2 (outer boundary halfedge)

    heh = _m.halfedge_handle(_eh, 
                             _m.is_boundary(_m.halfedge_handle(_eh,1)));
    
    assert( _m.is_boundary( _m.next_halfedge_handle( heh ) ) );
    assert( _m.is_boundary( _m.prev_halfedge_handle( heh ) ) );

    vh1 = _m.to_vertex_handle( _m.next_halfedge_handle( heh ) );
    vh2 = _m.to_vertex_handle( heh );
    vh3 = _m.from_vertex_handle( heh );
    vh4 = _m.from_vertex_handle( _m.prev_halfedge_handle( heh ));
    
    P1  = _m.point(vh1);
    P2  = _m.point(vh2);
    P3  = _m.point(vh3);
    P4  = _m.point(vh4);
    
    vhl = _m.add_vertex(zero);
    vhr = _m.add_vertex(zero);

    _m.property(vp_pos_, vhl ) = (P1 + real_t(16.0)*P2 + real_t(10.0)*P3) * _1over27;
    _m.property(vp_pos_, vhr ) = (real_t(10.0)*P2 + real_t(16.0)*P3 + P4) * _1over27;
    _m.property(ep_nv_, _eh).first  = vhl;
    _m.property(ep_nv_, _eh).second = vhr; 
  }


  void boundary_split( MeshType& _m, const typename MeshType::FaceHandle& _fh )
  {
    assert( _m.is_boundary(_fh) );

    typename MeshType::VertexHandle     vhl, vhr;
    typename MeshType::FaceEdgeIter     fe_it;
    typename MeshType::HalfedgeHandle   heh;

    // find boundary edge
    for( fe_it=_m.fe_iter( _fh ); fe_it && !_m.is_boundary( fe_it ); ++fe_it );

    // use precomputed, already inserted but not linked vertices
    vhl = _m.property(ep_nv_, fe_it).first;
    vhr = _m.property(ep_nv_, fe_it).second;

    /*
    //       *---------*---------*
    //      / \       / \       / \
    //     /   \     /   \     /   \
    //    /     \   /     \   /     \
    //   /       \ /       \ /       \
    //  *---------*--#---#--*---------*
    //                
    //  ^         ^  ^   ^  ^         ^
    //  P1        P2 pl  pr P3        P4
    */
    // get halfedge pointing from P2 to P3 (inner boundary halfedge)

    heh = _m.halfedge_handle(fe_it, 
                             _m.is_boundary(_m.halfedge_handle(fe_it,0)));

    typename MeshType::HalfedgeHandle pl_P3;

    // split P2->P3 (heh) in P2->pl (heh) and pl->P3
    boundary_split( _m, heh, vhl );         // split edge
    pl_P3 = _m.next_halfedge_handle( heh ); // store next halfedge handle
    boundary_split( _m, heh );              // split face

    // split pl->P3 in pl->pr and pr->P3
    boundary_split( _m, pl_P3, vhr );
    boundary_split( _m, pl_P3 );

    assert( _m.is_boundary( vhl ) && _m.halfedge_handle(vhl).is_valid() );
    assert( _m.is_boundary( vhr ) && _m.halfedge_handle(vhr).is_valid() );
  }

  void boundary_split(MeshType& _m, 
                      const typename MeshType::HalfedgeHandle& _heh,
                      const typename MeshType::VertexHandle& _vh)
  {
    assert( _m.is_boundary( _m.edge_handle(_heh) ) );

    typename MeshType::HalfedgeHandle 
      heh(_heh),
      opp_heh( _m.opposite_halfedge_handle(_heh) ),
      new_heh, opp_new_heh;
    typename MeshType::VertexHandle   to_vh(_m.to_vertex_handle(heh));
    typename MeshType::HalfedgeHandle t_heh;
    
    /*
     *            P5
     *             *
     *            /|\
     *           /   \
     *          /     \
     *         /       \
     *        /         \
     *       /_ heh  new \
     *      *-----\*-----\*\-----*
     *             ^      ^ t_heh
     *            _vh     to_vh
     *
     *     P1     P2     P3     P4
     */
    // Re-Setting Handles
    
    // find halfedge point from P4 to P3
    for(t_heh = heh; 
        _m.next_halfedge_handle(t_heh) != opp_heh; 
        t_heh = _m.opposite_halfedge_handle(_m.next_halfedge_handle(t_heh)))
    {}
    
    assert( _m.is_boundary( t_heh ) );

    new_heh     = _m.new_edge( _vh, to_vh );
    opp_new_heh = _m.opposite_halfedge_handle(new_heh);

    // update halfedge connectivity

    _m.set_next_halfedge_handle(t_heh,   opp_new_heh); // P4-P3 -> P3-P2
    // P2-P3 -> P3-P5
    _m.set_next_halfedge_handle(new_heh, _m.next_halfedge_handle(heh));    
    _m.set_next_halfedge_handle(heh,         new_heh); // P1-P2 -> P2-P3
    _m.set_next_halfedge_handle(opp_new_heh, opp_heh); // P3-P2 -> P2-P1

    // both opposite halfedges point to same face
    _m.set_face_handle(opp_new_heh, _m.face_handle(opp_heh));

    // let heh finally point to new inserted vertex
    _m.set_vertex_handle(heh, _vh); 

    // let heh and new_heh point to same face
    _m.set_face_handle(new_heh, _m.face_handle(heh));

    // let opp_new_heh be the new outgoing halfedge for to_vh 
    // (replaces for opp_heh)
    _m.set_halfedge_handle( to_vh, opp_new_heh );

    // let opp_heh be the outgoing halfedge for _vh
    _m.set_halfedge_handle( _vh, opp_heh );
  }

  void boundary_split( MeshType& _m, 
                       const typename MeshType::HalfedgeHandle& _heh)
  {
    assert( _m.is_boundary( _m.opposite_halfedge_handle( _heh ) ) );

    typename MeshType::HalfedgeHandle 
      heh(_heh),
      n_heh(_m.next_halfedge_handle(heh));

    typename MeshType::VertexHandle   
      to_vh(_m.to_vertex_handle(heh));

    typename MeshType::HalfedgeHandle 
      heh2(_m.new_edge(to_vh,
                       _m.to_vertex_handle(_m.next_halfedge_handle(n_heh)))),
      heh3(_m.opposite_halfedge_handle(heh2));

    typename MeshType::FaceHandle
      new_fh(_m.new_face()),
      fh(_m.face_handle(heh));
    
    // Relink (half)edges    

#define set_next_heh set_next_halfedge_handle
#define next_heh next_halfedge_handle

    _m.set_face_handle(heh,  new_fh);
    _m.set_face_handle(heh2, new_fh);
    _m.set_next_heh(heh2, _m.next_heh(_m.next_heh(n_heh)));
    _m.set_next_heh(heh,  heh2);
    _m.set_face_handle( _m.next_heh(heh2), new_fh);

    // _m.set_face_handle( _m.next_heh(_m.next_heh(heh2)), new_fh);

    _m.set_next_heh(heh3,                           n_heh);
    _m.set_next_heh(_m.next_halfedge_handle(n_heh), heh3);
    _m.set_face_handle(heh3,  fh);
    // _m.set_face_handle(n_heh, fh);

    _m.set_halfedge_handle(    fh, n_heh);
    _m.set_halfedge_handle(new_fh, heh);

#undef set_next_halfedge_handle
#undef next_halfedge_handle

  }

private:

  weights_t     weights_;
  OpenMesh::VPropHandleT< typename MeshType::Point >                    vp_pos_;
  OpenMesh::EPropHandleT< std::pair< typename MeshType::VertexHandle,
                                     typename MeshType::VertexHandle> > ep_nv_;
  OpenMesh::MPropHandleT< size_t >                                      mp_gen_;

  const real_t _1over3;
  const real_t _1over27;
};


//=============================================================================
} // END_NS_UNIFORM
} // END_NS_SUBDIVIDER
} // END_NS_OPENMESH
//=============================================================================
#endif // OPENMESH_SUBDIVIDER_UNIFORM_SQRT3T_HH
//=============================================================================