Python shapely.geometry.MultiLineString() Examples

def test_cover_geometry_empty_geoms(tiler):
    """Empty geometries should return empty iterators."""
    assert not cover_geometry(tiler, geometry.Point(), 0) == True
    assert not cover_geometry(tiler, geometry.Point(), [0, 1]) == True
    assert not cover_geometry(tiler, geometry.MultiPoint(), 0) == True
    assert not cover_geometry(tiler, geometry.MultiPoint(), [0, 1]) == True
    assert not cover_geometry(tiler, geometry.LineString(), 0) == True
    assert not cover_geometry(tiler, geometry.LineString(), [0, 1]) == True
    assert not cover_geometry(tiler, geometry.MultiLineString(), 0) == True
    assert not cover_geometry(tiler, geometry.MultiLineString(), [0, 1]) == True
    assert not cover_geometry(tiler, geometry.Polygon(), 0) == True
    assert not cover_geometry(tiler, geometry.Polygon(), [0, 1]) == True
    assert not cover_geometry(tiler, geometry.MultiPolygon(), 0) == True
    assert not cover_geometry(tiler, geometry.MultiPolygon(), [0, 1]) == True
    assert not cover_geometry(tiler, geometry.GeometryCollection(), 0) == True
    assert not cover_geometry(tiler, geometry.GeometryCollection(), [0, 1]) == True 

def _line_iterator(obj):
    if isinstance(obj, (sg.LineString)):
        yield obj
    elif isinstance(obj, (sg.MultiLineString)):
        for obj2 in obj.geoms:
            yield obj2
    elif isinstance(obj, (sg.Polygon)):
        yield sg.LineString(obj.exterior)
        for obj2 in obj.interiors:
            yield sg.LineString(obj2)
    elif isinstance(obj, (sg.MultiPolygon)):
        for obj2 in obj.geoms:
            yield sg.LineString(obj2.exterior)
            for obj3 in obj2.interiors:
                yield sg.LineString(obj3)
    else:
        try:
            tup = tuple(obj)
        except TypeError:
            raise TypeError('Could not use type %s' % type(obj))
        else:
            for obj2 in tup:
                for line in _line_iterator(obj2):
                    yield line 

def _round_multilinestring_coords(mls, precision):
    """
    Round the coordinates of a shapely MultiLineString to some decimal precision.

    Parameters
    ----------
    mls : shapely.geometry.MultiLineString
        the MultiLineString to round the coordinates of
    precision : int
        decimal precision to round coordinates to

    Returns
    -------
    shapely.geometry.MultiLineString
    """
    return MultiLineString([_round_linestring_coords(ls, precision) for ls in mls]) 

def alpha_shape(points, alpha):
    """
    Compute the alpha shape (concave hull) of a set
    of points.
    @param points: Iterable container of points.
    @param alpha: alpha value to influence the
        gooeyness of the border. Smaller numbers
        don't fall inward as much as larger numbers.
        Too large, and you lose everything!
    """
    if len(points) < 4:
        # When you have a triangle, there is no sense
        # in computing an alpha shape.
        return geometry.MultiPoint(list(points)).convex_hull

    coords = np.array([point.coords[0] for point in points])
    tri = Delaunay(coords)
    triangles = coords[tri.vertices]
    a = ((triangles[:,0,0] - triangles[:,1,0]) ** 2 + (triangles[:,0,1] - triangles[:,1,1]) ** 2) ** 0.5
    b = ((triangles[:,1,0] - triangles[:,2,0]) ** 2 + (triangles[:,1,1] - triangles[:,2,1]) ** 2) ** 0.5
    c = ((triangles[:,2,0] - triangles[:,0,0]) ** 2 + (triangles[:,2,1] - triangles[:,0,1]) ** 2) ** 0.5
    s = ( a + b + c ) / 2.0
    areas = (s*(s-a)*(s-b)*(s-c)) ** 0.5
    circums = a * b * c / (4.0 * areas)
    filtered = triangles[circums < (1.0 / alpha)]
    edge1 = filtered[:,(0,1)]
    edge2 = filtered[:,(1,2)]
    edge3 = filtered[:,(2,0)]
    edge_points = np.unique(np.concatenate((edge1,edge2,edge3)), axis = 0).tolist()
    m = geometry.MultiLineString(edge_points)
    triangles = list(polygonize(m))
    return cascaded_union(triangles), edge_points 

def test_cover_geometry_multilinestring2(tiler, mls):
    """A MultiLineString geometry."""
    tiles = [tile for tile in cover_geometry(tiler, mls, 12)]
    assert len(tiles) == 47
    assert set(tiles) == {(2036, 1420, 12), (2036, 1421, 12), (2036, 1422, 12), (2036, 1423, 12),
                          (2036, 1424, 12), (2037, 1421, 12), (2037, 1422, 12), (2038, 1420, 12),
                          (2038, 1421, 12), (2038, 1422, 12), (2038, 1423, 12), (2038, 1424, 12),
                          (2039, 1420, 12), (2039, 1421, 12), (2039, 1422, 12), (2039, 1423, 12),
                          (2039, 1424, 12), (2040, 1420, 12), (2040, 1422, 12), (2040, 1424, 12),
                          (2041, 1420, 12), (2041, 1422, 12), (2041, 1424, 12), (2042, 1420, 12),
                          (2042, 1421, 12), (2042, 1422, 12), (2042, 1423, 12), (2042, 1424, 12),
                          (2043, 1420, 12), (2044, 1420, 12), (2044, 1421, 12), (2044, 1422, 12),
                          (2044, 1423, 12), (2044, 1424, 12), (2045, 1420, 12), (2046, 1420, 12),
                          (2046, 1421, 12), (2046, 1422, 12), (2046, 1423, 12), (2047, 1420, 12),
                          (2047, 1423, 12), (2047, 1424, 12), (2048, 1420, 12), (2048, 1421, 12),
                          (2048, 1422, 12), (2048, 1423, 12), (2048, 1424, 12)} 

def shapely_to_paths(g):
    if isinstance(g, geometry.Point):
        return []
    elif isinstance(g, geometry.LineString):
        return [list(g.coords)]
    elif isinstance(g, (geometry.MultiPoint, geometry.MultiLineString, geometry.MultiPolygon, geometry.collection.GeometryCollection)):
        paths = []
        for x in g:
            paths.extend(shapely_to_paths(x))
        return paths
    elif isinstance(g, geometry.Polygon):
        paths = []
        paths.append(list(g.exterior.coords))
        for interior in g.interiors:
            paths.extend(shapely_to_paths(interior))
        return paths
    else:
        raise Exception('unhandled shapely geometry: %s' % type(g)) 

def shapely_to_paths(g):
    if isinstance(g, geometry.Point):
        return []
    elif isinstance(g, geometry.LineString):
        return [list(g.coords)]
    elif isinstance(g, (geometry.MultiPoint, geometry.MultiLineString, geometry.MultiPolygon, geometry.collection.GeometryCollection)):
        paths = []
        for x in g:
            paths.extend(shapely_to_paths(x))
        return paths
    elif isinstance(g, geometry.Polygon):
        paths = []
        paths.append(list(g.exterior.coords))
        for interior in g.interiors:
            paths.extend(shapely_to_paths(interior))
        return paths
    else:
        raise Exception('unhandled shapely geometry: %s' % type(g)) 

def _explode_lines(shape):
    """
    Return a list of LineStrings which make up the shape.
    """

    if shape.geom_type == 'LineString':
        return [shape]

    elif shape.geom_type == 'MultiLineString':
        return shape.geoms

    elif shape.geom_type == 'GeometryCollection':
        lines = []
        for geom in shape.geoms:
            lines.extend(_explode_lines(geom))
        return lines

    return [] 

def paths_to_shapely(paths):
    # TODO: Polygons for closed paths?
    return geometry.MultiLineString(paths) 

def paths_to_shapely(paths):
    return geometry.MultiLineString(paths) 

def polygon_splits(n, x, y, r, b):
    lines = []
    for i in range(n):
        t = 2 * pi / n * (i + 0.5)
        lines.append([(x, y), (x + r * cos(t), y + r * sin(t))])
    return geometry.MultiLineString(lines).buffer(b) 

def parse_gpx(self, filename):
        tracks_layer = fiona.open(filename, layer='tracks')
        feature      = tracks_layer[0]
        self.points  = feature['geometry']['coordinates']
        self.track   = self.check_track(MultiLineString(self.points))
        self.name    = feature['properties']['name'] 

def create(cls, geom, face_centers):
        """
        Create from existing facets and just select the ones that lie inside this polygon.
        """
        if isinstance(geom, (LineString, MultiLineString)):
            return HybridGeometry(geom, set())
        faces = tuple(
            set(np.argwhere(shapely_to_mpl(subgeom).contains_points(face_centers)).flatten())
            for subgeom in assert_multipolygon(geom)
        )
        return HybridGeometry(geom, tuple(f for f in faces if f)) 

def process(self, processors) -> MultiLineString:
        """
        Generate the compound geometries based on the provided processors. Sub-class must
        override this function in their implementation.
        :param processors: list of processors
        :return: compound geometries
        """ 

def __plot_edges(self, to_plot: List, edges: List, vertices: List, label: str):
        lines = []
        for edge_idx in to_plot:
            edge = edges[edge_idx]
            start_vertex = vertices[edge.start]
            end_vertex = vertices[edge.end]
            # Note it could be that the edge is supposed to be parabola (edge.is_linear
            # will be false), but in our case we always have boxes with 90 degree
            # corners. If it's a parabola then the focus is one of these corners, and by
            # drawing a line instead of a parabola we at worse cut through this point,
            # which is fine.
            lines.append(
                geometry.LineString(
                    [[start_vertex.X, start_vertex.Y], [end_vertex.X, end_vertex.Y]]
                )
            )
        merged_line = ops.linemerge(geometry.MultiLineString(lines))
        kwargs = {"label": label, "alpha": 0.5, "color": self.__colour_mapping[label]}
        # Merged line is either a MultiLineString which means we need to draw multiple
        # lines, or it is a LineString which means we only need to draw one.
        if isinstance(merged_line, geometry.MultiLineString):
            for line in merged_line:
                xs, ys = self.__simplify_outlines(line)
                self.__ax.plot(xs, ys, **kwargs)
                kwargs.pop(
                    "label", None
                )  # Only pass label once for single legend entry
        else:
            xs, ys = self.__simplify_outlines(merged_line)
            self.__ax.plot(xs, ys, **kwargs) 

def generate() -> MultiLineString:
    """
    A python module called by the `script` command must implement this function, which takes
    no arguments and return a MultiLineString object.
    :return: resulting MultiLineString
    """
    segs = []

    # horizontal segments
    for i in range(math.floor(N / 2)):
        for j in range(M):
            append_maybe([(i, j), (i + 1, j)], segs)

    # add half horizontal segments
    if N % 2 == 0:
        for j in range(M):
            append_maybe([((N / 2) - 1, j), (N / 2, j)], segs)

    # add vertical segments
    for i in range(math.ceil(N / 2)):
        for j in range(M - 1):
            append_maybe([(i, j), (i, j + 1)], segs)

    half_mls = MultiLineString(segs)
    other_mls = affinity.translate(affinity.scale(half_mls, -1, 1, origin=(0, 0)), N - 1, 0)

    return ops.linemerge(ops.unary_union([half_mls, other_mls])) 

def as_mls(self) -> MultiLineString:
        return MultiLineString([as_vector(line) for line in self.lines]) 

def test_geometry_coords_rounding():
    # test the rounding of geometry coordinates
    precision = 3

    shape1 = Point(1.123456, 2.123456)
    shape2 = ox.utils_geo.round_geometry_coords(shape1, precision)

    shape1 = MultiPoint([(1.123456, 2.123456), (3.123456, 4.123456)])
    shape2 = ox.utils_geo.round_geometry_coords(shape1, precision)

    shape1 = LineString([(1.123456, 2.123456), (3.123456, 4.123456)])
    shape2 = ox.utils_geo.round_geometry_coords(shape1, precision)

    shape1 = MultiLineString(
        [
            [(1.123456, 2.123456), (3.123456, 4.123456)],
            [(11.123456, 12.123456), (13.123456, 14.123456)],
        ]
    )
    shape2 = ox.utils_geo.round_geometry_coords(shape1, precision)

    shape1 = Polygon([(1.123456, 2.123456), (3.123456, 4.123456), (6.123456, 5.123456)])
    shape2 = ox.utils_geo.round_geometry_coords(shape1, precision)

    shape1 = MultiPolygon(
        [
            Polygon([(1.123456, 2.123456), (3.123456, 4.123456), (6.123456, 5.123456)]),
            Polygon([(16.123456, 15.123456), (13.123456, 14.123456), (12.123456, 11.123456)]),
        ]
    )
    shape2 = ox.utils_geo.round_geometry_coords(shape1, precision) 

def round_geometry_coords(shape, precision):
    """
    Round the coordinates of a shapely geometry to some decimal precision.

    Parameters
    ----------
    shape : shapely.geometry.geometry
        the geometry to round the coordinates of; must be one of {Point,
        MultiPoint, LineString, MultiLineString, Polygon, MultiPolygon}
    precision : int
        decimal precision to round coordinates to

    Returns
    -------
    shapely.geometry.geometry
    """
    if isinstance(shape, Point):
        return _round_point_coords(shape, precision)

    elif isinstance(shape, MultiPoint):
        return _round_multipoint_coords(shape, precision)

    elif isinstance(shape, LineString):
        return _round_linestring_coords(shape, precision)

    elif isinstance(shape, MultiLineString):
        return _round_multilinestring_coords(shape, precision)

    elif isinstance(shape, Polygon):
        return _round_polygon_coords(shape, precision)

    elif isinstance(shape, MultiPolygon):
        return _round_multipolygon_coords(shape, precision)

    else:
        raise TypeError(f"cannot round coordinates of unhandled geometry type: {type(shape)}") 

def _split_multigeom_row(gdf_row, geom_col):
    new_rows = []
    if isinstance(gdf_row[geom_col], MultiPolygon) \
            or isinstance(gdf_row[geom_col], MultiLineString):
        new_polys = _split_multigeom(gdf_row[geom_col])
        for poly in new_polys:
            row_w_poly = gdf_row.copy()
            row_w_poly[geom_col] = poly
            new_rows.append(row_w_poly)
    return pd.DataFrame(new_rows) 

def redistribute_vertices(geom, dist):
    """
    Redistribute the vertices on a projected LineString or MultiLineString.

    The distance argument is only approximate since the total distance of the
    linestring may not be a multiple of the preferred distance. This function
    works on only (Multi)LineString geometry types.

    Parameters
    ----------
    geom : shapely.geometry.LineString or shapely.geometry.MultiLineString
        a Shapely geometry (should be projected)
    dist : float
        spacing length along edges. Units are same as the geom: degrees for
        unprojected geometries and meters for projected geometries. The
        smaller the dist value, the more points are created.

    Returns
    -------
    list or shapely.geometry.MultiLineString
        the redistributed vertices as a list if geom is a LineString or
        MultiLineString if geom is a MultiLineString
    """
    if geom.geom_type == "LineString":
        num_vert = int(round(geom.length / dist))
        if num_vert == 0:
            num_vert = 1
        return [geom.interpolate(float(n) / num_vert, normalized=True) for n in range(num_vert + 1)]
    elif geom.geom_type == "MultiLineString":
        parts = [redistribute_vertices(part, dist) for part in geom]
        return type(geom)([p for p in parts if not p])
    else:
        raise ValueError(f"unhandled geometry {geom.geom_type}") 

def _lines_only(shape):
    """
    Extract the lines (LineString, MultiLineString) from any geometry. We
    expect the input to be mostly lines, such as the result of an intersection
    between a line and a polygon. The main idea is to remove points, and any
    other geometry which might throw a wrench in the works.
    """

    lines = _explode_lines(shape)
    if len(lines) == 1:
        return lines[0]
    else:
        return MultiLineString(lines) 

def assert_multilinestring(geometry: Union[LineString, MultiLineString, GeometryCollection]) -> List[LineString]:
    """
    given a LineString or MultiLineString, return a list of LineStrings
    :param geometry: a LineString or a MultiLineString
    :return: a list of LineStrings
    """
    if geometry.is_empty:
        return []
    if isinstance(geometry, LineString):
        return [geometry]
    return [geom for geom in geometry.geoms if isinstance(geom, LineString)] 

def create_full(cls, geom, vertices_offset, faces_offset):
        """
        Create by triangulating a polygon and adding the resulting facets to the total list.
        """
        if isinstance(geom, (LineString, MultiLineString, Point)):
            return HybridGeometry(geom, set()), np.empty((0, 2), dtype=np.int32), np.empty((0, 3), dtype=np.uint32)

        vertices = deque()
        faces = deque()
        faces_i = deque()
        for subgeom in assert_multipolygon(geom):
            new_vertices, new_faces = triangulate_polygon(subgeom)
            new_faces += vertices_offset
            vertices.append(new_vertices)
            faces.append(new_faces)
            faces_i.append(set(range(faces_offset, faces_offset+new_faces.shape[0])))
            vertices_offset += new_vertices.shape[0]
            faces_offset += new_faces.shape[0]

        if not vertices:
            return HybridGeometry(geom, set()), np.empty((0, 2), dtype=np.int32), np.empty((0, 3), dtype=np.uint32)

        vertices = np.vstack(vertices)
        faces = np.vstack(faces)

        return HybridGeometry(geom, tuple(faces_i)), vertices, faces 

def to_shapely(self):
        return geometry.MultiLineString(self.paths) 

def GetUsCanadaBorder():
  """Gets the US/Canada border as a |shapely.MultiLineString|."""
  global _uscanada_border
  if _uscanada_border is None:
    kml_file = os.path.join(CONFIG.GetFccDir(), USCANADA_BORDER_FILE)
    lines = _ReadKmlBorder(kml_file)
    _uscanada_border = ops.unary_union(lines.values())
  return _uscanada_border 

def _ReadKmlBorder(kml_path, root_id='Placemark'):
  """Gets the border defined in a KML.

  Args:
    kml_path: The path name to the border file KML or KMZ.
    root_id_zone: The root id defininig a zone. Usually it is 'Placemark'.

  Returns:
    A dictionary of |shapely| LineString keyed by their names.
  """
  if kml_path.endswith('kmz'):
    with zipfile.ZipFile(kml_path) as kmz:
      kml_name = [info.filename for info in kmz.infolist()
                  if os.path.splitext(info.filename)[1] == '.kml'][0]
      with kmz.open(kml_name) as kml_file:
        root = parser.parse(kml_file).getroot()
  else:
    with open(kml_path, 'r') as kml_file:
      root = parser.parse(kml_file).getroot()

  tag = root.tag[:root.tag.rfind('}') + 1]
  linetrings_dict = {}
  for element in root.findall('.//' + tag + root_id):
    # Ignore nested root_id within root_id
    if element.find('.//' + tag + root_id) is not None:
      continue
    name = element.name.text
    linestrings = [
        _GetLineString(l)
        for l in element.findall('.//' + tag + 'LineString')
    ]
    if not linestrings:
      continue
    if len(linestrings) == 1:
      linestring = linestrings[0]
    else:
      linestring = sgeo.MultiLineString(linestrings)

    linetrings_dict[name] = linestring

  return linetrings_dict 

def test_cover_geometry_multilinestring3(tiler, mls):
    """A MultiLineString geometry."""
    tiles = [tile for tile in cover_geometry(tiler, mls, [7, 8])]
    assert len(tiles) == 2
    assert set(tiles) == {(63, 44, 7), (64, 44, 7)} 

def test_cover_geometry_multilinestring1(tiler, mls):
    """A MultiLineString geometry."""
    tiles = [tile for tile in cover_geometry(tiler, mls, 8)]
    assert len(tiles) == 4
    assert set(tiles) == {(127, 88, 8), (127, 89, 8), (128, 88, 8), (128, 89, 8)} 

def mls():
    return geometry.shape({'type': 'MultiLineString', 'coordinates': (((-0.9709167480468749, 35.18840002173177), (-0.9832763671875, 34.83296838321102), (-0.954437255859375, 35.02774729487063), (-0.850067138671875, 35.02662273458687), (-0.833587646484375, 34.838604318635014), (-0.8074951171874999, 35.185032937998294)), ((-0.0823974609375, 35.191766965947394), (-0.111236572265625, 34.83071390101431), (0.048065185546875, 34.829586636768205), (0.067291259765625, 35.192889249680945), (-0.06591796875, 35.191766965947394)), ((-0.28564453125, 35.19064466671118), (-0.34606933593749994, 34.84085858477277), (-0.18402099609375, 34.831841149828676)), ((-0.48614501953124994, 35.183910545750834), (-0.517730712890625, 34.854382885097905), (-0.391387939453125, 34.84874803007872)), ((-0.553436279296875, 35.183910545750834), (-0.7347106933593749, 35.185032937998294), (-0.751190185546875, 35.05922870088872), (-0.591888427734375, 35.0502352513963), (-0.758056640625, 35.03336986422378), (-0.767669677734375, 34.86227103378598), (-0.597381591796875, 34.85550980979319)))})