I have a 2D model of a road that consists of separate rectangles whose ends are approximately lined up. Below an example of a situation where the road turns slightly.

Ends of two roads approximately lined up

There are thousands of rectangles, and that's a problem for later processing. I need to join ends of most of the roads such that they form a connected mesh, like so:

Ends (edges) of two roads joined

I can't simply perform a Remove Doubles on all the geometry, because this is done from a script, and the risk of joining vertices in a way that doesn't result in sensible geometry is too high. For example below only the bottom vertices of the edge have been joined. In intersections where roads overlap, completely unrelated roads might become connected in broken ways.

Ends of two roads where just one vertex of the edge has been joined

Thus what I need is an operation that joins each edge pair whose ends are sufficiently close to each other. Is there any way easier than writing that entirely by myself? If I were to write that myself, what would be the best strategy?

The end goal is to create a model that can be 3D printed nicely. I extrude the 2D roads, resulting in a tactile map that you get from Touch Mapper. Thousands of objects is a bit much for slicers, which basically need to first separate by loose parts, and then perform an union between every resulting object.

Below is an example of roads in a real-word city centre. Pedestrian roads (being edited in the screenshot) are extruded higher than other roads. I'd like to join eg. the edge whose vertices are hilighted.

Screenshot of geometry in a busy intersection

Edit: final code

Below is my final code. Changes to Chebhou's answer are:

  • Speed up edge proximity searching with mathutils.kdtree.KDTree. Important with over a few thousand edges.
  • Merge all verts in one operation using bmesh.ops.weld_verts(). With 1000 verts to join, this shrinks execution time from several seconds to almost nothing.
  • Filter out weird edges, and precompute some values in filter_edges().
  • Don't merge edges in places where there more than 2 edges meet. This is a good heuristic that avoids all roads becoming a single object, which would not be desirable when 3D printing, because objects shouldn't intersect themselves.

Mesh is assumed to be in object mode initially.

def join_matching_edges(ob):
  lt = 0.1  # length difference + -
  dt = 0.1  # max distance
  at = 0.5  # max sin(angle)  (30°)

  bpy.context.scene.objects.active = ob
  bpy.ops.object.mode_set(mode = 'EDIT')
  from math import sin
  bm = bmesh.from_edit_mesh( bpy.context.object.data )

  center   = lambda e : ( e.verts[0].co + e.verts[1].co ) / 2
  length   = lambda e : ( e.verts[0].co - e.verts[1].co ).length
  dist     = lambda v1, v2: (  v2 -  v1 ).length
  sinAngle = lambda e1, e2: abs(sin((e1.verts[1].co - e1.verts[0].co).angle(e2.verts[1].co - e2.verts[0].co)))

  def point_between_edge_neighbor_verts(e):
      # Return middle of the verts adjacent to the edge
      verts = []
      for v in e.verts:
          for linked_e in v.link_edges:
              verts.extend((vv for vv in linked_e.verts if vv != e.verts[0] and vv != e.verts[1]))
      if len(verts) != 2:
          #print("edge has non-2 adjacent verts: " + str(len(verts)))
          return None
      return ((verts[0].co[0] + verts[1].co[0]) / 2, \
              (verts[0].co[1] + verts[1].co[1]) / 2, \
              (verts[0].co[2] + verts[1].co[2]) / 2)

  class CEdge:
      def __init__(self, e, center, into_edge):
          self.e = e
          self.center = center(e)
          self.length = length(e)
          self.into_edge = into_edge
          self.welded = False

  def filter_edges(edges):
      out = []
      for e in edges:
          if len(e.link_faces) != 1:
          point_between_edges = point_between_edge_neighbor_verts(e)
          if not point_between_edges:
          vector_into_edge_face = center(e) - mathutils.Vector(point_between_edges)
          out.append(CEdge(e, center, vector_into_edge_face / vector_into_edge_face.length))
      return out
  candidate_edges = filter_edges(bm.edges)

  # Index edges into search tree
  kd = mathutils.kdtree.KDTree(len(candidate_edges))
  for i, ce in enumerate(candidate_edges):
      kd.insert(ce.center, i)

  to_weld = {}
  for i, ce in enumerate(candidate_edges[:-1]):
      if ce.welded:
      ce.welded = True
      lmin = ce.length - lt
      lmax = ce.length + lt
      matches = []
      for (_co, oe_index, _dist) in kd.find_range(ce.center, dt):
          oe = candidate_edges[oe_index]
          if not oe.welded and lmin < oe.length < lmax and sinAngle(ce.e, oe.e) < at:
              turn_angle = ce.into_edge.angle(-oe.into_edge)
              if turn_angle > math.pi * 0.6: # pi * 0.5 is 90%
                  #print("not merging edges (%s, %s) pointing to opposite directions, angle is %f" % (ce.e, oe.e, turn_angle))
              oe.welded = True

      if len(matches) == 1:
          # Join nothing where >2 ways meet, else all roads in the scene may become joined and cross itself
          ev1, ev2 = ce.e.verts[:]
          oev1, oev2 = matches[0].e.verts[:]
          if dist(ev1.co, oev1.co) < dist(ev1.co, oev2.co) :
              if ev1 != oev1 : to_weld[ev1] = oev1
              if ev2 != oev2 : to_weld[ev2] = oev2
          else :
              if ev1 != oev2 : to_weld[ev1] = oev2
              if ev2 != oev1 : to_weld[ev2] = oev1

  print("%s: melding %d out of %d edges" % (ob.name, len(to_weld) / 2, len(bm.edges)))
  bmesh.ops.weld_verts(bm, targetmap = to_weld)
  bmesh.update_edit_mesh(bpy.context.object.data ,True)
  bpy.ops.object.mode_set(mode = 'OBJECT')
  • $\begingroup$ You could try select one of these edges and select others by direction from Select Similar operator (Shift+G). Adjusting Threshold in the operator settings, you could select only desired edges. Then remove doubles. This way won't really work if there're important edges similarly directed; probably you'll need to control visibility of the objects. $\endgroup$ – Mr Zak Apr 8 '16 at 17:26
  • $\begingroup$ There will be all sorts of edges in there, but Select Similar could still be useful as part of the solution. Thanks. $\endgroup$ – Samuli Pahaoja Apr 8 '16 at 17:52

This is modified version of TLousky answer in this script we use three factors to decide merging edges :

  • distance between their centers.
  • difference in length.
  • the angle between these edges.

you can drop any of them if it's not necessary in your scene, and this is the script :

import bpy, bmesh
from math import sin

bm = bmesh.from_edit_mesh( bpy.context.object.data )

lt = 0.2 #length difference + -
dt = 0.2  #max distance 
at = 0.16  #max sin(angle)  (30°)

center  = lambda e : ( e.verts[0].co + e.verts[1].co ) / 2
length  = lambda e : ( e.verts[0].co - e.verts[1].co ).length
dist    = lambda v1, v2: (  v2 -  v1 ).length
sinAngle   = lambda e1, e2: abs(sin((e1.verts[1].co - e1.verts[0].co).angle(e2.verts[1].co - e2.verts[0].co))) 

i =-1
while i < len(bm.edges)-1:
    i +=1
    e = bm.edges[i]
    i = e.index 
    l = length(e)
    lmin = l-lt
    lmax = l+lt
    c = center(e)
    for oe in  bm.edges[i+1:]:
        if  lmin < length(oe) < lmax  and dist( center(oe), c) < dt and sinAngle(e, oe)<at :
            ev1, ev2 = e.verts[:]
            oev1, oev2 = oe.verts[:]
            print( ev2, oev2 )
            if dist(ev1.co, oev1.co ) < dist(ev1.co, oev2.co ) :
                if ev1 != oev1 : bmesh.ops.pointmerge( bm, verts = [ ev1, oev1 ], merge_co = ev1.co )
                if ev2 != oev2 : bmesh.ops.pointmerge( bm, verts = [ ev2, oev2 ], merge_co = ev2.co )
            else :
                if ev1 != oev2 : bmesh.ops.pointmerge( bm, verts = [ ev1, oev2 ], merge_co = ev1.co )
                if ev2 != oev1 : bmesh.ops.pointmerge( bm, verts = [ ev2, oev1 ], merge_co = ev2.co )

bmesh.update_edit_mesh(bpy.context.object.data ,True)
  • 1
    $\begingroup$ Works perfectly, after I added bm.edges.ensure_lookup_table() after bm's initialization. I added that into the answer's code. Thank you! $\endgroup$ – Samuli Pahaoja Apr 26 '16 at 5:36

This script will delete edges that are closer to each other than the threshold specified edgeThreshold if both their vertices are closer than the vertThreshold.

I use two thresholds in order to avoid iterating over all the edges for each edge, and filtering faraway edges (will make a big difference if the mesh is large and has many edges).

It assumes your mesh is in edit mode.

enter image description here

import bpy, bmesh

bm = bmesh.from_edit_mesh( bpy.context.object.data )

mat           = bpy.context.object.matrix_world
vertThreshold = 0.3
edgeThreshold = 0.4

calcEdgeCenter = lambda e, mat: ( mat * e.verts[0].co + mat * e.verts[1].co ) / 2
dist           = lambda v1, v2, mat: ( mat * v2 - mat * v1 ).length

removedEdges = []
for i, e in enumerate( bm.edges ):
    if i in removedEdges: continue

    # Find current edge's center point
    curEdgeCenter = calcEdgeCenter( e, mat )

    # Find close edges (whose center is closer than the edgeThreshold)
    otherEdges = [ 
        oe for oe in bm.edges if oe.index != e.index and \
        oe.index not in removedEdges and \
        dist( curEdgeCenter, calcEdgeCenter( oe, mat ), mat ) < edgeThreshold

    eV1, eV2 = e.verts[:]

    for oe in otherEdges:
        oeIdx = oe.index
        oeV1, oeV2 = oe.verts[:]

        # An edge will be considered redundant if both its verts are closer to 
        # the current edge's verts. In that case merge the edges' vertices
        if dist( eV1.co, oeV1.co, mat ) < vertThreshold and dist( eV2.co, oeV2.co, mat ) < vertThreshold:
            if eV1.index != oeV1.index:
                bmesh.ops.pointmerge( bm, verts = [ eV1, oeV1 ], merge_co = mat * eV1.co )
            if eV2.index != oeV2.index:
                bmesh.ops.pointmerge( bm, verts = [ eV2, oeV2 ], merge_co = mat * eV2.co )
            removedEdges += [ oeIdx, i ]
        elif dist( eV1.co, oeV2.co, mat ) < vertThreshold and dist( eV2.co, oeV1.co, mat ) < vertThreshold:
            if eV1.index != oeV2.index:
                bmesh.ops.pointmerge( bm, verts = [ eV1, oeV2 ], merge_co = mat * eV1.co )
            if eV2.index != oeV1.index:
                bmesh.ops.pointmerge( bm, verts = [ eV2, oeV1 ], merge_co = mat * eV2.co )
            removedEdges += [ oeIdx, i ]
  • $\begingroup$ If you read my question carefully, I mention that simply removing double verts isn't reliable enough. The geometry is unpredictable, and completely unrelated verts may happen to be arbitrary close to each other. $\endgroup$ – Samuli Pahaoja Apr 8 '16 at 17:42
  • $\begingroup$ Missed that @SamuliKärkkäinen. Will update my answer to devise a better trick. $\endgroup$ – TLousky Apr 8 '16 at 17:45
  • $\begingroup$ @SamuliKärkkäinen, updated my answer with a new script that should do the trick. $\endgroup$ – TLousky Apr 8 '16 at 18:38
  • 1
    $\begingroup$ Outstanding answer! For me O(n*n) is too slow, but I'll see about speeding that up using KDTree. $\endgroup$ – Samuli Pahaoja Apr 9 '16 at 7:01
  • $\begingroup$ To make the script work on my scene, I needed to remove "mat" from the "merge_co = mat * eVx.co"'s. After that things work with low thresholds, but when I increase the thresholds, I get every time a crash, or "BMesh data of type BMVert has been removed" on "if dist( eV1.co, oeV1.co, mat ) < vertThreshold and dist( eV2.co, oeV2.co, mat ) < vertThreshold:". In the beginning of the inner loop oe.is_valid is sometimes false, but skipping those oe's doesn't end crashing. I can't figure this out. $\endgroup$ – Samuli Pahaoja Apr 9 '16 at 19:24

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.