add new tool refs #16051

tools/shapes/poly2edgedata.py

@@ -0,0 +1,283 @@
+#!/usr/bin/env python
+# Eclipse SUMO, Simulation of Urban MObility; see https://eclipse.dev/sumo
+# Copyright (C) 2010-2025 German Aerospace Center (DLR) and others.
+# This program and the accompanying materials are made available under the
+# terms of the Eclipse Public License 2.0 which is available at
+# https://www.eclipse.org/legal/epl-2.0/
+# This Source Code may also be made available under the following Secondary
+# Licenses when the conditions for such availability set forth in the Eclipse
+# Public License 2.0 are satisfied: GNU General Public License, version 2
+# or later which is available at
+# https://www.gnu.org/licenses/old-licenses/gpl-2.0-standalone.html
+# SPDX-License-Identifier: EPL-2.0 OR GPL-2.0-or-later
+
+# @file    poly2edgedata.py
+# @author  Jakob Erdmann
+# @date    2025-01-23
+
+"""
+Transform polygons with params into edgedata with attributes
+For each polygon a unique edge is selected that gives the best geometrical match
+
+The following syntax is supported in the patch file (one patch per line):
+ 
+# lines starting with '#' are ignored as comments
+# rev overrides the reverse edge of EDGEID to be REVEDGEID
+rev EDGEID REVEDGEID
+# edg overrides the edge to assign for POLYID to be EDGEID
+edg POLYID EDGEID
+# dat overrrides the data attribute ATTR for POLYID to take on value VALUE
+dat POLYID ATTR VALUE
+
+any ID or VALUE may bet set to 'None' to signify that
+- a reverse edge should not be assigned
+- a polygon should not be mapped
+- data should be ignored
+"""
+
+from __future__ import print_function
+from __future__ import absolute_import
+import os
+import sys
+from math import fabs, degrees
+from collections import defaultdict
+THIS_DIR = os.path.dirname(os.path.abspath(__file__))
+sys.path.append(os.path.join(THIS_DIR, '..'))
+import sumolib  # noqa
+import sumolib.geomhelper as gh
+
+PATCH_NONE = 'None'
+
+def get_options(args=None):
+    op = sumolib.options.ArgumentParser(description="Transform polygons with params into edgedata with attributes")
+
+    op.add_option("polyfiles", nargs='+', category="input", type=op.file_list,
+                  help="List of polygon files to convert")
+    op.add_option("-n", "--netfile", category="input", required=True, type=op.net_file,
+                  help="Network file on which to map the polygons")
+    op.add_option("-p", "--patchfile", category="input",  type=op.file,
+                    help="Load a file with patches to apply during import")
+    op.add_option("-r", "--radius", type=float, default=20,
+                  help="radius for finding edges near polygons")
+    op.add_option("--shapecut", type=float, default=40,
+                  help="Shorten polygon and edge shapes to FLOAT to increase robustness of angular comparison")
+    op.add_option("--angle-tolerance", type=float, default=20, dest="atol",
+                  help="Match polygons and edges if their angle differs by no more than DEGRESS")
+    op.add_option("--min-length", type=float, default=10, dest="minLength",
+                  help="minimum edge length that may be mapped to")
+    op.add_option("-s", "--split-attributes", dest="splitAttrs",
+                  help="If a reverse edge is found, split the values of the given attribute list among edge and reverse edge")
+    op.add_option("-S", "--nosplit-attributes", dest="noSplitAttrs",
+                  help="If a reverse edge is found, split the values of all attributes except the given attribute list among edge and reverse edge")  # noqa
+    op.add_option("-f", "--filter", dest="filter", 
+                  help="Read a list of triplets ATTR,MIN,MAX and only keep polygons where value ATTR is within [MIN,MAX]")
+    op.add_option("-b", "--begin", default=0, type=op.time,
+                  help="edgedata interval begin time")
+    op.add_option("-e", "--end", default="1:0:0:0", type=op.time,
+                  help="edgedata interval end time)")
+    op.add_option("-o", "--output-file", category="output", dest="outfile", required=True, type=op.file,
+                  help="output file")
+
+    try:
+        options = op.parse_args()
+    except (NotImplementedError, ValueError) as e:
+        print(e, file=sys.stderr)
+        sys.exit(1)
+
+    options.splitAttrs = set(options.splitAttrs.split(',')) if options.splitAttrs else []
+    options.noSplitAttrs = set(options.noSplitAttrs.split(',')) if options.noSplitAttrs else []
+    tuples = options.filter.split(',') if options.filter else []
+    options.filter = {} # attr -> (min, max)
+    for i in range(0, len(tuples), 3):
+        options.filter[tuples[i]] = (float(tuples[i + 1]), float(tuples[i + 2]))
+    return options
+
+
+def hasReverse(edge):
+    for cand in edge.getToNode().getOutgoing():
+        if cand.getToNode() == edge.getFromNode():
+            return True
+    return False
+
+
+def readPatches(net, pfile):
+    patchEdg = {}  # polyID->edge
+    patchRev = {}  # forwardEdge->reverseEdge
+    patchDat = defaultdict(lambda: {}) # polyID->attr->data (-1 ignores)
+    with open(pfile) as pf:
+        for line in pf:
+            items = line.split()
+            patchtype = items[0]
+            if patchtype == "rev":
+                edgeID, reverseID = items[1:]
+                patchRev[edgeID] = reverseID
+                test = net.getEdge(edgeID)
+            elif patchtype == "edg":
+                polyID, edgeID = items[1:]
+                patchEdg[polyID] = edgeID
+            elif patchtype == "dat":
+                polyID, attrName, value = items[1:]
+                patchDat[polyID][attrName] = value
+            elif patchtype == "#":
+                # comment 
+                continue
+            else:
+                print("unknown patchtype '%s'" % patchtype)
+    return patchEdg, patchRev, patchDat
+
+
+def main(options):
+    net = sumolib.net.readNet(options.netfile)
+    patchEdg, patchRev, patchDat = readPatches(net, options.patchfile)
+    usedEdges = set() # do not assign different polygons/counts to the same edge
+    scut = options.shapecut
+
+    with open(options.outfile, 'w') as foutobj:
+        foutobj.write('<meandata>\n')
+        foutobj.write('    <interval begin="%s" end="%s" id="%s">\n' % (
+            options.begin, options.end, options.polyfiles[0]))
+        for fname in options.polyfiles:
+            for poly in sumolib.xml.parse(fname, 'poly'):
+
+                shape = []
+                for lonlat in poly.shape.split():
+                    lon,lat = lonlat.split(',')
+                    shape.append(net.convertLonLat2XY(float(lon), float(lat)))
+                shapelen = gh.polyLength(shape)
+                cx, cy = gh.positionAtShapeOffset(shape, shapelen / 2)
+                edges = net.getNeighboringEdges(cx, cy, options.radius)
+                if not edges:
+                    print("No edges near %.2f,%.2f (poly %s)" % (cx, cy, poly.id), file=sys.stderr)
+                    continue
+                edges = [(e, d) for e, d in edges if e.allows("passenger")]
+                if not edges:
+                    print("No car edges near %.2f,%.2f (poly %s)" % (cx, cy, poly.id), file=sys.stderr)
+                    continue
+                edges = [(e, d) for e, d in edges if e.getLength() >= options.minLength]
+                if not edges:
+                    print("No long edges near %.2f,%.2f (poly %s)" % (cx, cy, poly.id), file=sys.stderr)
+                    continue
+
+                if shapelen < scut:
+                    polyAngle = gh.angleTo2D(shape[0], shape[-1])
+                else:
+                    polyAngle = gh.angleTo2D(gh.positionAtShapeOffset(shape, shapelen / 2 - scut / 2),
+                                             gh.positionAtShapeOffset(shape, shapelen / 2 + scut / 2))
+
+                cands = []
+                for e, d in edges:
+                    if e.getLength() < scut:
+                        angle = gh.angleTo2D(e.getFromNode().getCoord(), e.getToNode().getCoord())
+                        revAngle = gh.angleTo2D(e.getToNode().getCoord(), e.getFromNode().getCoord())
+                    else:
+                        eShape = e.getShape()
+                        offset = gh.polygonOffsetWithMinimumDistanceToPoint((cx, cy), eShape)
+                        offset1 = max(0, offset - scut / 2)
+                        offset2 = min(gh.polyLength(eShape), offset + scut / 2)
+                        angle = gh.angleTo2D(gh.positionAtShapeOffset(eShape, offset1),
+                                             gh.positionAtShapeOffset(eShape, offset2))
+                        revAngle = gh.angleTo2D(gh.positionAtShapeOffset(eShape, offset2),
+                                             gh.positionAtShapeOffset(eShape, offset1))
+                    if degrees(fabs(polyAngle - angle)) < options.atol or degrees(fabs(polyAngle - revAngle)) < options.atol:
+                        cands.append(e)
+                edges = cands
+                if not edges:
+                    print("No edges with angle %.2f found near %.2f,%.2f (poly %s)" % (
+                        degrees(polyAngle), cx, cy, poly.id), file=sys.stderr)
+                    continue
+
+                bestDist = 1e10
+                bestEdge = None
+                bestReverse = None
+                for e in edges:
+                    if shapelen > e.getLength():
+                        maxDist = max([gh.distancePointToPolygon(xy, shape) for xy in e.getShape()])
+                    else:
+                        maxDist = max([gh.distancePointToPolygon(xy, e.getShape()) for xy in shape])
+                    #if poly.id == "241517.0":
+                    #    print(e.getID(), maxDist, "shapelen", shapelen, "eLen", e.getLength())
+                    if maxDist < bestDist:
+                        bestDist = maxDist
+                        bestEdge = e
+                # apply edge patch
+                if poly.id in patchEdg:
+                    if patchEdg[poly.id] == PATCH_NONE:
+                        continue
+                    bestEdge = net.getEdge(patchEdg[poly.id])
+
+                if bestEdge in usedEdges:
+                    patchInfo = " (was patched)" if poly.id in patchEdg else ""
+                    print("Duplicate assignment to edge %s from poly %s%s" % (
+                        bestEdge.getID(), poly.id, patchInfo),  file=sys.stderr)
+                    continue
+                # find opposite direction for undivided road
+                for e in edges:
+                    if (e.getFromNode() == bestEdge.getToNode() and
+                        e.getToNode() == bestEdge.getFromNode()):
+                        bestReverse = e
+                        break
+                # apply revers edge patch
+                if bestEdge.getID() in patchRev:
+                    revID = patchRev.get(bestEdge.getID())
+                    if revID != PATCH_NONE:
+                        bestReverse = net.getEdge(revID)
+                    else:
+                        bestReverse = None
+                # find opposite direction for divided road
+                elif bestReverse is None:
+                    bestDist = 1e10
+                    bestAngle = gh.angleTo2D(bestEdge.getFromNode().getCoord(), bestEdge.getToNode().getCoord())
+                    for e in edges:
+                        if hasReverse(e):
+                            continue
+                        reverseAngle = gh.angleTo2D(e.getToNode().getCoord(), e.getFromNode().getCoord())
+                        if degrees(fabs(bestAngle - reverseAngle)) < 20:
+                            maxDist = max([gh.distancePointToPolygon(xy, bestEdge.getShape()) for xy in e.getShape()])
+                            if maxDist < bestDist:
+                                bestDist = maxDist
+                                bestReverse = e
+
+                if bestReverse in usedEdges:
+                    patchInfo = " (was patched)" if bestEdge.getID() in patchRev else ""
+                    print("Duplicate assignment to reverse edge %s from poly %s%s" % (
+                        bestReverse.getID(), poly.id, patchInfo),  file=sys.stderr)
+                    continue
+
+                attrs = 'polyID="%s"' % poly.id
+                skip = False
+                for param in poly.param:
+                    value = param.value
+                    if poly.id in patchDat:
+                        if param.key in patchDat[poly.id]:
+                            value = patchDat[poly.id][param.key]
+                            #print("patched %s to %s" % (param.value, value))
+                            if value == PATCH_NONE:
+                                continue
+                    if param.key in options.filter:
+                        if (float(value) < options.filter[param.key][0] or
+                            float(value) > options.filter[param.key][1]):
+                            skip = True
+                    if bestReverse is not None and (param.key in options.splitAttrs
+                            or options.noSplitAttrs and param.key not in options.noSplitAttrs):
+                        try:
+                            value = float(value) / 2
+                        except:
+                            pass
+                    attrs += ' %s="%s"' % (param.key, value)
+
+                if skip:
+                    continue
+                comment = ''
+                if bestReverse is not None:
+                    usedEdges.add(bestReverse)
+                    comment = '  <!-- reverse: %s -->' % bestEdge.getID()
+                    foutobj.write('        <edge id="%s" %s/> %s\n' % (bestReverse.getID(), attrs, comment))
+                    comment = '  <!-- reverse: %s -->' % bestReverse.getID()
+                usedEdges.add(bestEdge)
+                foutobj.write('        <edge id="%s" %s/> %s\n' % (bestEdge.getID(), attrs, comment))
+            foutobj.write('    </interval>\n')
+            foutobj.write('</meandata>\n')
+
+
+if __name__ == "__main__":
+    main(get_options())