Source code for ddg.datastructures.halfedge.modify

"""Utils to modify a given surface or its cells"""
import logging
from itertools import cycle, islice

import numpy as np

import ddg
from ddg import nonexact
from ddg.datastructures.halfedge.get import (
    count_edges_in_loop,
    get_edge_loop,
    get_edges,
    get_vertex,
    get_vertices,
    in_edges,
    is_boundary_edge,
    is_boundary_vertex,
    out_edges,
    single_edges,
)

logger = logging.getLogger(__name__)

################################
# joining, attaching and glueing
################################


[docs]def fill_hole(e): """ Fills a hole, i.e. a loop of boundary edges, with a face. Parameters ---------- e : edge An interior halfedge belonging to the loop of boundary edges that is to fill. Returns ------- face Raises ------ ValueError If there already is a face and no hole to fill. """ if not is_boundary_edge(e): raise ValueError( "No hole to fill. Argument e already has a face. " + "An edge of a boundary loop was expected." ) s = e.surf # face = s.add_polygon(count_edges_in_loop(e)) # face_edge = face.edge # edge_loop = list(get_edge_loop(e)) # for edge in edge_loop: # s.glue(edge.opp, face_edge) # face_edge = face_edge.nex f = s.faces() f.edge = e for edge in get_edge_loop(e): edge.face = f return f
[docs]def join_neighbouring_faces(e): """ Joins two adjacent faces by removing the edge they have in common. Parameters ---------- e : edge The edge to be removed of its corresponding surface. Returns ------- None """ if e.opp.face: remove_face(e.opp.face) # assign e.face to all edges and remove references pointing to e or e.opp for edge in get_edge_loop(e.opp): edge.face = e.face if e.face.edge == e: e.face.edge = e.nex if e.head.edge == e: e.head.edge = e.nex.opp if e.tail.edge == e.opp: e.tail.edge = e.opp.nex.opp # set new edge succesor and predecessor in_head = e.opp.pre out_head = e.nex in_tail = e.pre out_tail = e.opp.nex _set_next(in_head, out_head) _set_next(in_tail, out_tail) # remove edges e.surf.edges.remove(e.opp) e.surf.edges.remove(e)
[docs]def attach_pyramid(e): """ Attaches a pyramid to the loop of boundary edges along a given one. One new triangle each is glued to every edge of the edge loop of consecutive boundary edges along `e` and as well glued to the previously added triangle, such that they all share the same tip and form a pyramid. Parameters ---------- e : edge An edge of the boundary loop along which a pyramid will be attached. Returns ------- vertex The newly added vertex which is the peak of the attached pyramid. Raises ------ ValueError If the edge does not lie opposite to a boundary loop. """ if not is_boundary_edge(e): raise ValueError( "The edge does not lie opposite to a boundary loop " "and thus no pyramid can be attached." ) s = e.surf boundary_loop = list(get_edge_loop(e)) for i, (edge) in enumerate(boundary_loop): triangle = s.add_face(3) s.glue(edge.opp, triangle.edge) # face_edge = face_edge.nex if i == 0: first_vertical_edge = triangle.edge.pre else: s.glue(vertical_edge, triangle.edge.pre) vertical_edge = triangle.edge.nex if i == len(boundary_loop) - 1: peak = first_vertical_edge.tail s.glue(first_vertical_edge, vertical_edge) return peak
[docs]def bridge_loops(e1, e2): """ Creates a bridging loop between two not yet connected loops of boundary edges or components of a surface by gluing a circular strip of new edges and quad-faces in. The loops of edges that shall be bridged should have the same number of edges. The bridging starts with a quad connecting 'edge1' and 'edge2'. Parameters ---------- e1 : edge An interior halfedge opposite to a boundary loop that shall be bridged. e2 : edge An interior halfedge opposite to the boundary loop that shall be bridged to the first. Returns ------- None Raises ------ ValueError If the loops don't count the same number of edges. """ s = e1.surf s2 = e2.surf if s is not s2: raise ValueError( "The edges do not belong to the same surface and thus cannot be bridged." ) if not is_boundary_edge(e1.opp): raise ValueError( "The first edge parameter does not lie opposite to a " "boundary loop and thus cannot be bridged." ) if not is_boundary_edge(e2.opp): raise ValueError( "The second edge parameter does not lie opposite to a" " boundary loop and thus cannot be bridged." ) if count_edges_in_loop(e1.opp) is not count_edges_in_loop(e2.opp): raise ValueError( "The edge loops do not count the same number of edges" " and thus cannot be bridged. \n " "loop1: {} edges, loop2: {} edges".format( count_edges_in_loop(e1.opp), count_edges_in_loop(e2.opp) ) ) top_edge = e1 bottom_edge = e2 top_edges = [top_edge] bottom_edges = [bottom_edge] top_edge = top_edge.opp.pre.opp bottom_edge = bottom_edge.opp.nex.opp while top_edge is not e1 and bottom_edge is not e2: top_edges.append(top_edge) bottom_edges.append(bottom_edge) top_edge = top_edge.opp.pre.opp bottom_edge = bottom_edge.opp.nex.opp for i, (top_edge, bottom_edge) in enumerate(zip(top_edges, bottom_edges)): quad = s.add_face(4) s.glue(top_edge, quad.edge) s.glue(bottom_edge, quad.edge.nex.nex) if i == 0: left_vertical_edge = quad.edge.nex else: s.glue(vertical_edge, quad.edge.nex) vertical_edge = quad.edge.pre if i == len(top_edges) - 1: s.glue(left_vertical_edge, vertical_edge)
[docs]def extrude(f, remove=True): """ Extrude a face of a surface, adding a new face between each halfedge. Parameters ---------- f : face A face of a surface. remove : bool (default=True) If True, remove the face f after extrusion. Returns ------- The extruded face. """ face_edge = f.edge surf = f.surf if remove: remove_face(f) if is_boundary_edge(face_edge): start_edge = face_edge else: start_edge = face_edge.opp n = count_edges_in_loop(start_edge) extrude_face = surf.add_face(n) bridge_loops(extrude_face.edge, start_edge.opp) return extrude_face
[docs]def zip_digon(f): r""" Zips a digon to an edge Parameters ---------- f : face A digon Returns ------- f.edge.opp : edge Raises ------ ValueError if parameter f is not a digon """ e = f.edge surf = e.surf if e.nex == e.pre: # digon e_nex = e.nex e_opp = e.opp e.head.edge = e_nex.opp e.tail.edge = e_opp e_opp.opp = e_nex.opp e_nex.opp.opp = e_opp surf.edges.remove(e) surf.edges.remove(e_nex) surf.faces.remove(f) return e_opp else: raise ValueError("Face f is not a digon")
[docs]def join_coplanar_faces(surf, co_attr="co", atol=None, rtol=None): """ Investigates whether a surface has triangulated, planar faces and if needed the triangulation will be removed, i.e. triangles will be joined to a face with more than three vertices. Parameters ---------- surf : surface The halfedge surface. co_attr : str The attribute name that contains the coordinate attribute of the vertices. atol, rtol : float or None (default=None) Tolerances used to determine whether two faces are coplanar. This function uses the global tolerance defaults if `atol` or `rtol` are set to None. See :py:mod:`ddg.nonexact` for details. """ if not all((hasattr(v, str(co_attr)) for v in surf.verts)): raise ValueError("None or not all vertices have a coordinate attribute") coordinates = getattr(surf.verts, co_attr) def is_almost_coplanar(edge): A = [ np.array(coordinates[edge.pre.head]) - np.array(coordinates[edge.head]), np.array(coordinates[edge.head]) - np.array(coordinates[edge.nex.head]), np.array(coordinates[edge.head]) - np.array(coordinates[edge.opp.pre.pre.head]), ] return nonexact.isclose(np.linalg.det(A), 0.0, atol=atol, rtol=rtol) for e in [e for e in surf.edges if is_almost_coplanar(e)]: if e in surf.edges: join_neighbouring_faces(e)
def _insert_edge(v1, v2): """ Inserts a new edge between given vertices of a surface. Parameters ---------- v1, v2: vertex The vertices that will be connected by a new edge. They need to belong to the same surface. Returns ------- new_e : edge The halfedge of the newly inserted edge that points from v1 towards v2, i.e. with `new_e.tail = v1` and `new_e.head = v2`. Notes ----- Attention: This method only inserts the two opposing halfedges between the given vertices. It does not link the halfedges to any other edge of the surface, thus the references `pre`, `nex` and `face` are None and need to be set in the further process. Raises ------ ValueError If the given vertices do not belong to the same surface. """ surf = v1.surf surf2 = v2.surf if surf is not surf2: raise ValueError( "The vertices do not belong to the same surface " "and thus cannot be connected by an edge." ) # insert two new halfedges to form an edge between vertex v1 and v2 new_e = surf.edges() new_e_opp = surf.edges() # and establish references new_e.opp = new_e_opp new_e_opp.opp = new_e new_e.head = v2 new_e_opp.head = v1 if v2.edge is None: v2.edge = new_e if v1.edge is None: v1.edge = new_e.opp return new_e def _link_edge_in(e_pre, e, e_nex): """ Links an edge between a given previous and next halfedge by establishing all references between them. Parameters ---------- e_pre : edge The halfedge that shall become `e.pre`. e : edge The halfedge to be linked inbetween `e_pre` and `e_nex`. e_nex : edge The halfedge that shall become `e.nex`. Returns ------- None Notes ----- By linking `e` to given previous and next halfedges, all references are updated accordingly (e.g. `e_pre.nex` will be set to `e`). Attention: BUT the reference of former partners of `e`, `e_pre` and `e_nex` to them will be dissolved, meaning `e.nex.pre`, `e.pre.nex`, `e_nex.pre.nex` and `e_pre.nex.pre` will be set to None! This may cause Errors if not taken care of! """ _set_next(e, e_nex) _set_previous(e, e_pre) def _set_next(e, nex): """ Sets the next halfedge of an edge to be the given one. Parameters ---------- e : edge The edge whose 'e.nex' shall become 'nex'. nex : edge The edge that shall become the next edge of 'e'. Returns ------- None Notes ----- By linking `e` to the given next halfedge all references are updated accordingly (e.g. `nex.pre` will be set to `e`). Attention: BUT the reference of former partners of `e` and `nex` to them will be dissolved, meaning `e.nex.pre`, and `nex.pre.nex` will be set to None! This may cause Errors if not taken care of! """ if e.nex is nex: return if e.nex is not None: e.nex.pre = None if nex is not None: if nex.pre is not None: nex.pre.nex = None nex.pre = e e.nex = nex def _set_previous(e, pre): """ Sets the previous halfedge of an edge to be the given one. Parameters ---------- e : edge The edge whose 'e.pre' shall become 'pre'. pre : edge The edge that shall become the previous edge of 'e'. Returns ------- None Notes ----- By linking `e` to the given previous halfedge all references are updated accordingly (e.g. `pre.nex` will be set to `e`). Attention: BUT the reference of former partners of `e` and `pre` to them will be dissolved, meaning `e.pre.nex` and `pre.nex.pre` will be set to None! This may cause Errors if not taken care of! """ if e.pre is pre: return if e.pre is not None: e.pre.nex = None if pre is not None: if pre.nex is not None: pre.nex.pre = None pre.nex = e e.pre = pre #################################### # separating, splitting and removing ####################################
[docs]def remove_face(f): """ Removes a specified face from its surface. The corresponding halfedges will be retained and their face set to be None. Parameters ---------- f : face The face to be removed of its corresponding surface. Returns ------- None """ for e in get_edges(f): e.face = None s = f.surf s.faces.remove(f)
[docs]def remove_vertex(v): """Remove a vertex from a surface. Also removes all incoming and outgoing edges and any faces that these edges belong to. Parameters ---------- v : vertex The vertex to be removed from the surface it belongs to. """ surf = v.surf star = list(in_edges(v)) for e in star: remove_edge(e) surf.verts.remove(v)
[docs]def remove_edge(e): """Remove an edge from a surface. Also removes faces that `e` or `e.opp` belong to. Parameters ---------- e : edge """ surf = e.surf if e.face: remove_face(e.face) if e.opp.face: remove_face(e.opp.face) e.pre.nex = e.opp.nex e.opp.nex.pre = e.pre e.nex.pre = e.opp.pre e.opp.pre.nex = e.nex e.tail.edge = e.pre e.head.edge = e.opp.pre surf.edges.remove(e.opp) surf.edges.remove(e)
[docs]def split_face_at(f, v1, v2): """ Insert a new edge between two vertices and thereby split a face into two. Parameters ---------- f : face The face that will be split into two faces. v1, v2 : vertex The vertices between a new edge will be inserted to split the face. They need to be vertices in the boundary of `f`. Returns ------- new_edge : edge Returns the just added halfedge which is pointing towards v2, i.e. with `edge.head == v2` and `edge.tail == v1`. """ # find the incoming and outgoing edges of v1 and v2 in the specified face f logger.debug("face to split:" + str(f.index)) logger.debug("its vertices:" + str([v.index for v in get_vertices(f)])) logger.debug("v1,v2 for split:" + str([v1.index, v2.index])) tmp_in = [] tmp_out = [] for v in [v1, v2]: for edge in in_edges(v): if edge.face is f: tmp_in.append(edge) if edge.opp.face is f: tmp_out.append(edge.opp) in_v1, in_v2 = tmp_in out_v1, out_v2 = tmp_out assert in_v1 is not None and out_v1 is not None assert in_v2 is not None and out_v2 is not None if out_v1 is in_v2 or out_v2 is in_v1: raise ValueError( "No degenerate split allowed. The given vertices are " "neighbouring, thus the face cannot be split." ) # insert new edge and its opposite between v1 and v2 s = f.surf new_edge = _insert_edge(v1, v2) _link_edge_in(in_v1, new_edge, out_v2) _link_edge_in(in_v2, new_edge.opp, out_v1) # the new edge will belong to the old face f new_edge.face = f f.edge = new_edge # insert second face and set references to belong to new_edge.opp # and all edges of its edgeloop new_face = s.faces() new_face.edge = new_edge.opp new_edge.opp.face = new_face for e in get_edge_loop(new_edge.opp): e.face = new_face return new_edge
[docs]def stellar_subdivide(f): """ Adds a new vertex and subdivides the given face by stellar subdivision. The given face `f` subdivides into number of boundary edges many triangular faces. The former face is removed. Parameters ---------- f : face Returns ------- peak : vertex The newly added vertex. See Also -------- datastructures.halfedge.modify.subdivide : Subdivides a half-edge surface's faces by bisecting the edges and connecting the new vertices. """ edge = f.edge remove_face(f) peak = attach_pyramid(edge) return peak
[docs]def subdivide_edge(e, number_of_segments): """ Subdivides an edge into the given number of segments. Parameters ---------- e : edge A halfedge of the edge that will be subdivided. number_of_segments : int The number of segments the edge will be subdivided into. Returns ------- list containing * segment of newly subdivided edge based at `e.tail` * segment of newly subdivided edge based at `e.opp.tail` * List of all vertices that were added by subdivision. """ if number_of_segments == 1: return [e, e.opp, []] surf = e.surf head = e.head first_segment = e first_opp_segment = e.opp edge = e added_vertices = [] for _ in range(number_of_segments - 1): # insert a new vertex and set its references v = surf.verts() added_vertices.append(v) v.surf = surf v.edge = edge # insert an edge between v and e.head and set face references new_e = _insert_edge(v, head) if edge.face is not None: new_e.face = edge.face if edge.opp.face is not None: new_e.opp.face = edge.opp.face # update references of involved old and new vertices and edges if head.edge is edge: head.edge = new_e edge.head = v _link_edge_in(edge, new_e, edge.nex) _link_edge_in(edge.opp.pre, new_e.opp, edge.opp) edge = new_e # store first of the opposite segments if _ is number_of_segments - 2: first_opp_segment = new_e.opp # first_opp_segment = first_segment # for _ in range(number_of_segments-1): # first_opp_segment = first_opp_segment.nex # first_opp_segment = first_opp_segment.opp return [first_segment, first_opp_segment, added_vertices]
[docs]def contract_edge(e): r"""Contract an edge and remove the edge and its opposite from the surface. Note that this can alter the manifold property if an adjacent face is a triangle or `e.head` and `e.tail` are boundary vertices:: *-----* * /| | /| * |e | -> *--* | \| | \| *-----* * *-----*-----* * * | | | |\ /| | |e | -> | * | | | | |/ \| *-----*-----* * * Parameters ---------- e : edge to be contracted Returns ------- edge, edge e.pre and e.opp.pre of the contracted edges. Be careful: these edges might also have been deleted if they coincided with e or e.opp """ surf = e.surf e_opp = e.opp e_pre, e_nex = e.pre, e.nex e_opp_pre, e_opp_nex = e.opp.pre, e.opp.nex if e.nex == e: if e.face: surf.faces.remove(e.face) e.head.edge = None else: e_pre.nex = e_nex e_nex.pre = e_pre e.head.edge = e_pre if e.face: e.face.edge = e_pre if e_opp_nex == e_opp: if e_opp.face: surf.faces.remove(e_opp.face) e_opp.head.edge = None else: e_opp_pre.nex = e_opp_nex e_opp_nex.pre = e_opp_pre e_opp.head.edge = e_opp_pre if e_opp.face: e_opp.face.edge = e_opp_pre if e.head != e_opp.head: for in_edge in in_edges(e.head): in_edge.head = e_opp.head surf.verts.remove(e.head) surf.edges.remove(e) surf.edges.remove(e_opp) return e_pre, e_opp_pre
[docs]def contract_face(f): r"""Contract a face to a vertex. Note that this can alter the manifold property if neighboring faces of `f` are triangles or opposite sides of `f` are part of a boundary:: *-----* /| | * | f | -> *--* \| | *-----* *-----*-----*-----* * * | | | | |\ /| | | f | | -> | * | | | | | |/ \| *-----*-----*-----* * * Parameters ---------- f : face to be contracted Returns ------- vertex the face was contracted to """ fe = f.edge v = fe.tail remove_face(f) for e in list(get_edge_loop(fe)): contract_edge(e) return v
[docs]def subdivide(surface, steps=1, co_attr="co"): """ Subdivides a half-edge surface's faces by bisecting the edges and connecting the new vertices. Bisects every edge of the surface into two segments and a new vertex. Then the vertices of any face form a cycle in which the old and new vertices alternate. The new vertices are connected in this order. This procedure is repeated steps many times. The number of new faces is exponential in steps. Parameters ---------- surf : ddg.halfege.Surface The half-edge surface. steps: int (default=1) Number of repeated subdivisions. co_attr : str or None (default='co') A string can be given representing the attribute name that contains the coordinate attribute of the vertices. Then the vertex subdividing the edge is set to its midpoint. If None is given the subdivision is purely combinatorial. Returns ------- None See Also -------- ddg.halfedge.modify.stellar_subdivide : Adds a new vertex and subdivides the given face by stellar subdivision. """ def basic_bisect(edge): _, _, x = subdivide_edge(edge, 2) x[0].old = False return x[0] if co_attr is None: bisect = basic_bisect else: def bisect(edge): # Calculate coordinates before edge is split and destroyed. co = (getattr(edge.tail, co_attr) + getattr(edge.head, co_attr)) / 2 x = basic_bisect(edge) setattr(x, co_attr, co) return x for _ in range(steps): surface.verts.add_attribute("old", default=True) # We add new faces in the loop, so take the list to ensure we only loop # over old faces. for face in list(surface.faces): # We filter out v := edge.head if v.old because v will be the tail # of the next edge in get_edge_loop(face.edge). By construction, # new and old faces will alternate. new_face_verts = list( map( lambda edge: bisect(edge) if edge.tail.old else edge.tail, filter( lambda edge: edge.head.old, # Similarly to before, take the list to make # sure we don't iterate over and mutate # face.edge at the same time. list(get_edge_loop(face.edge)), ), ) ) for v, w in zip(new_face_verts, islice(cycle(new_face_verts), 1, None)): split_face_at(face, v, w) del surface.verts.old
##################### # reverse_orientation #####################
[docs]def reverse_orientation(s): """ Function to reverse the orientation of edges of a surface. Parameters ---------- s: ddg.halfedge.Surface Surface to reverse the orientation on """ for e in single_edges(s): ehead = e.head etail = e.tail enex = e.nex epre = e.pre eoppnex = e.opp.nex eopppre = e.opp.pre e.nex = epre e.pre = enex e.head = etail e.opp.nex = eopppre e.opp.pre = eoppnex e.opp.head = ehead for v in s.verts: v.edge = v.edge.opp
[docs]def diagonals(hds, boundary_face_incidence=4, append_single_edges=True): r"""Creates two hds' with the diagonal combinatorics to the given hds. The given hds will be bi-colored and if a vertex of one color is adjacent to more than two vertices of the other color this vertex will be replaced by a face build from the surrounding vertices. In this sense two new halfedge data structures are build. The vertices inherit all (manually set) attributes of the given vertices an additionally have an attribute 'old_index' storing the previous vertex index. Parameters ---------- hds: ddg.datastructures.halfedge.Surface Surface to create the hds' from. It must be able to apply halfedge.set.bicolor_vertices to this. boundary_face_incidence: int (default=4) If a boundary vertex is incident to at least this many vertices of the other color, then a face will be added in the above sense. If it is adjacent to less, only the edges will be added, excluding the edge along the boundary of the original surface. append_single_edges: bool (default=True) If set to 'True' edges that do not belong to a face will be added (i.e., 'tentacles' at the corners of rectangular grid-like hds). This means the resulting hds might not be a manifold. Raises ------ ddg.datastructures.indexedfaceset.utils.BoundaryException If in the resulting surface two faces are only connected to each other by one vertex, like this:: /\/\ \/\/ Returns ------- hds1, hds2 The resulting hds' of diagonal combinatoric. """ verts_attr_list = [ attr for attr in get_vertex(hds).__dir__() if (attr[0] != "_" and attr not in ["edge", "surf", "index"]) ] verts_attr_dict = { v.index: {attr: getattr(v, attr) for attr in verts_attr_list} for v in hds.verts } temporary_color_attr = "temporary_color_attr" ddg.halfedge.set.bicolor_vertices(hds, color_attr=temporary_color_attr) diagonals = [None, None] faces = [[], []] edges = [set(), set()] # creating the face lists and edge sets for v in hds.verts: face = [e.head for e in out_edges(v)] k = len(face) # usual inner face or face at the boundary with enough vertices to build a face if ( not is_boundary_vertex(v) and k > 2 or is_boundary_vertex(v) and k >= boundary_face_incidence ): faces[getattr(v, temporary_color_attr) - 1].append([v.index for v in face]) # edges might not belong to a face so we remember them else: if k == 2: edges[getattr(v, temporary_color_attr) - 1].add( tuple(sorted((face[0].index, face[1].index))) ) else: for i in range(0, k): # exclude boundary edge if not is_boundary_vertex(face[i]) or not is_boundary_vertex( face[(i + 1) % k] ): edges[getattr(v, temporary_color_attr) - 1].add( tuple(sorted((face[i].index, face[(i + 1) % k].index))) ) # build hds' from the face lists and append single_edges if stated for i in range(2): if faces[i]: face_set = ddg.indexedfaceset.ifs.GeneralizedIndexedFaceSet(faces[i]) hds_i = ddg.indexedfaceset.utils.indexed_face_set_to_surface(face_set) else: hds_i = ddg.datastructures.halfedge.Surface() for attr in verts_attr_list: hds_i.verts.add_attribute(attr) hds_i.verts.add_attribute("old_index") if append_single_edges: # Tentacles are all the edges that are not covered by the edges of the faces if faces[i]: single_edges = edges[i] - face_set.edge_set() else: single_edges = edges[i] single_edge_vertices = { index: True for index in set([ind for tuple in edges[i] for ind in tuple]) } for v in hds_i.verts: for attr, value in verts_attr_dict[v.ifs_index].items(): setattr(v, attr, value) v.old_index = v.ifs_index if append_single_edges: if single_edge_vertices.get(v.ifs_index, False): single_edge_vertices[v.ifs_index] = v if append_single_edges: for edge in single_edges: # attach_tentacle requires new vertex as last second vertex keyword offset = -1 if single_edge_vertices[edge[0]] is True else 0 for j in range(2): if single_edge_vertices[edge[j]] is True: v = hds_i.verts() for attr, value in verts_attr_dict[edge[j]].items(): setattr(v, attr, value) v.old_index = edge[j] single_edge_vertices[edge[j]] = v hds_i.add_edge( single_edge_vertices[edge[0 + offset]], single_edge_vertices[edge[1 + offset]], ) if faces[i]: delattr(hds_i.verts, "ifs_index") diagonals[i] = hds_i delattr(hds.verts, temporary_color_attr) return diagonals[0], diagonals[1]