3

I'm currently working on a project involving creating an isochrone map of historical water travel. At the moment, I have a geodataframe containing a starting point on the river, a multipoint of the upstream and downstream endpoints that can be reached within the given timeframe, and a buffer around the starting point representing the distance that can be walked from the starting point in the given time. I also have a linestring of the river. I'm trying to create a polygon representing the total area that can be covered through a combination of boats and walking, so I basically need to create a buffer around the river that is widest at the starting point (the width of my current circular buffer) and tapers towards the ends. In the picture below, I have the river in blue and the points shown in black, and I need to generate a shape similar to the purple dotted outline. Drawing of a river with start and endpoints, as well as desired buffer

From previous research, I know I could create a buffer around the river, but only of uniform thickness. I could also generate a convex hull, but that would not conform well to the shape of the river. To make matters worse, I need to do this for about 400 separate locations, so I need a solution that can be coded. I'm fairly new to this, so apologies for anything I've overlooked!

2

1 Answer 1

0

Create a point every meter along the line from each end to the start, and buffer each point starting at one meter and finish at the middle points buffer width.

Then union all the buffers.

import numpy as np
from shapely import wkb
from shapely.ops import unary_union
import geopandas as gpd
import matplotlib.pyplot as plt

riverwkb = b'\x01\x02\x00\x00\x00H\x00\x00\x00\x10\xdf\xd0\xdf\xd99"A\xe8@%\x97!U\\Ac\x8aM\xba\xd89"A\xc4iN\xbe$U\\A\x18\xddZP\xdd9"Au\xbd\x16\x9c\'U\\A/\xd7\xf8\xa1\xe79"A\x91\xa6.U*U\\A\xec\'\x90\xa8\xef9"A\xd8d\xa7W-U\\A\x0b \xb8j\xfd9"A^\xe3\x0e\xec/U\\A2\x16j\x9d\x0e:"A""e\x122U\\AY\x0c\x1c\xd0\x1f:"AR\xf6\n\x144U\\A\xd3\xadJ\xdd/:"A\xac\x9f\x11_6U\\ALOy\xea?:"A\x07I\x18\xaa8U\\AdI\x17<J:"A\x8d\xc7\x7f>;U\\A\x19\x9c$\xd2N:"A>\x1bH\x1c>U\\A\x11\x9e\x9aaK:"A\xefn\x10\xfa@U\\A\xaf\xf6\t\xa6E:"A\xa0\xc2\xd8\xd7CU\\A\xfa\xa3\xfc\x0fA:"A\xe7\x80Q\xdaFU\\A\xf2\xa5r\x9f=:"A\x98\xd4\x19\xb8IU\\A\xf2\xa5r\x9f=:"AI(\xe2\x95LU\\A\x9f\xfa\xf5\xc4>:"A\xfa{\xaasOU\\AdI\x17<J:"A\x80\xfa\x11\x08RU\\A\x8a?\xc9n[:"AE9h.TU\\A\xb93\x05\x12p:"A\x1dcL\x9dUU\\A\xf0%\xcb%\x88:"A\nx\xbeTVU\\Ay\xc3\r\x14\x9f:"A\xf6\x8c0\x0cWU\\A\x02aP\x02\xb6:"AM7\xf2\x9eWU\\A8S\x16\x16\xce:"A\xa3\xe1\xb31XU\\AoE\xdc)\xe6:"A\xcf\xb6\x14{XU\\A\xa57\xa2=\xfe:"AQ67WYU\\A\x89~\xebv\x17;"A\xd3\xb5Y3ZU\\A\x12\x1c.e.;"AT5|\x0f[U\\A9\x12\xe0\x97?;"A\x19t\xd25]U\\AA\x10j\x08C;"A`2K8`U\\AA\x10j\x08C;"A=[t_cU\\A|\xc1H\x917;"AXD\x8c\x18fU\\AU\xcb\x96^&;"A\x87\x182\x1ahU\\A\'\xd7Z\xbb\x11;"A\xf6\xac\xc6\xadiU\\A\xf0\xe4\x94\xa7\xf9:"A\xf9\xab\x0bfkU\\A\xc2\xf0X\x04\xe5:"A\xd2\xd5\xef\xd4lU\\A\xe6\xa7\x99;\xcf:"A\x15\x95#\x1fnU\\A]\nWM\xb8:"A\x83)\xb8\xb2oU\\A.\x16\x1b\xaa\xa3:"A\xc6\xe8\xeb\xfcpU\\A\xff!\xdf\x06\x8f:"A4}\x80\x90rU\\A\xd0-\xa3cz:"A8|\xc5HtU\\A\xfc\xe2m\x0bh:"A\xd2\xe5\xba%vU\\A\xe5\xe8\xcf\xb9]:"A\xed\xce\xd2\xdexU\\A8\x94L\x94\\:"A\x9e"\x9b\xbc{U\\A8\x94L\x94\\:"APvc\x9a~U\\AO\x8e\xea\xe5f:"A\xd5\xf4\xca.\x81U\\A\xc11\x8f\x82s:"A[s2\xc3\x83U\\A\x95|\xc4\xda\x85:"A\x8aG\xd8\xc4\x85U\\A\xc4p\x00~\x9a:"A\xcd\x06\x0c\x0f\x87U\\AM\x0eCl\xb1:"AO\x86.\xeb\x87U\\A\xd6\xab\x85Z\xc8:"A\x10\xc6?Y\x88U\\A\r\x9eKn\xe0:"Agp\x01\xec\x88U\\A\x96;\x8e\\\xf7:"A\x92Eb5\x89U\\A\x1f\xd9\xd0J\x0e;"A\xe9\xef#\xc8\x89U\\A\xa8v\x139%;"A~Z\xd4\xec\x89U\\A1\x14V\'<;"A\xaa/56\x8aU\\A\xbb\xb1\x98\x15S;"A?\x9a\xe5Z\x8aU\\A\xf1\xa3^)k;"A\xd5\x04\x96\x7f\x8aU\\A\xd5\xea\xa7b\x84;"A\x96D\xa7\xed\x8aU\\A\x0b\xddmv\x9c;"AW\x84\xb8[\x8bU\\A\x94z\xb0d\xb3;"A\x82Y\x19\xa5\x8bU\\A\x1e\x18\xf3R\xca;"Aon\x8b\\\x8cU\\A\xa7\xb55A\xe1;"A[\x83\xfd\x13\x8dU\\A8Q\x02\xa0\xfb;"A\xdd\x02 \xf0\x8dU\\AnC\xc8\xb3\x13<"A\xc9\x17\x92\xa7\x8eU\\A\xf7\xe0\n\xa2*<"AK\x97\xb4\x83\x8fU\\A.\xd3\xd0\xb5B<"A\xf8\xeb7\xa9\x90U\\A\xb7p\x13\xa4Y<"AzkZ\x85\x91U\\A@\x0eV\x92p<"A\x92U-\x86\x92U\\Ao\x02\x925\x85<"Aj\x7f\x11\xf5\x93U\\A\xe8\xa3\xc0B\x95<"A/\xbeg\x1b\x96U\\A'
points = [b'\x01\x04\x00\x00\x00\x01\x00\x00\x00\x01\x01\x00\x00\x00\xf0\xe4\x94\xa7\xf9:"A\xf9\xab\x0bfkU\\A',
          b'\x01\x04\x00\x00\x00\x02\x00\x00\x00\x01\x01\x00\x00\x00\x9f\xfa\xf5\xc4>:"A\xfa{\xaasOU\\A\x01\x01\x00\x00\x008Q\x02\xa0\xfb;"A\xdd\x02 \xf0\x8dU\\A']

river = wkb.loads(riverwkb)
riverdf = gpd.GeoDataFrame(geometry=[river], crs="epsg:3006")

pointgeoms = [wkb.loads(p) for p in points]
pointdf = gpd.GeoDataFrame(data={"type":["start", "end"]}, geometry=pointgeoms, crs="epsg:3006")
#     type                                           geometry
# 0  start                MULTIPOINT (597372.827 7427501.594)
# 1    end  MULTIPOINT (597279.385 7427389.807, 597501.813...

#Multipoints to singlepoints
pointdf["points"] = pointdf.apply(lambda x: [pnt for pnt in x.geometry.geoms], axis=1)
pointdf = pointdf.explode(column="points", ignore_index=True).set_geometry("points").drop(columns="geometry").rename_geometry("geometry")


centerbuffer = pointdf.loc[pointdf["type"]=="start"].buffer(70)

#Calculate distance along the river for each of the three points
pointdf["distance_along_line"] = pointdf.apply(lambda x: river.line_locate_point(x.geometry), axis=1)
#     type                        geometry  distance_along_line
# 0  start  POINT (597372.827 7427501.594)           423.724419
# 1    end  POINT (597279.385 7427389.807)           203.194798
# 2    end  POINT (597501.813 7427639.752)           770.110404

d1 = abs(pointdf.distance_along_line.iloc[0]-pointdf.distance_along_line.iloc[1]) #Number of points between top point and middle
d2 = abs(pointdf.distance_along_line.iloc[0]-pointdf.distance_along_line.iloc[2]) #Number of points between top point and bottom

#Create a point every 1 m from top to center point
p1 = [river.line_interpolate_point(x) for x in np.arange(203, 423+1, 1)]
#Buffer them starting at 1 m up to 70 m
buffers1 = [pnt.buffer(dist) for pnt, dist in zip(p1, np.arange(start=1, stop=70, step=(70-1)/len(p1)))]

p2 = [river.line_interpolate_point(x) for x in np.arange(423, 770+1, 1)]
buffers2 = [pnt.buffer(dist) for pnt, dist in zip(p2, np.arange(start=70, stop=1, step=-(70-1)/len(p2)))]

uu = unary_union(buffers1+buffers2)
buffdf = gpd.GeoDataFrame(geometry=[uu], crs="epsg:3006")

#Plot
fig, ax = plt.subplots(figsize=(15, 15))

#buffdf.plot(ax=ax, color="green", alpha=0.2)
buffdf.boundary.plot(ax=ax, color="purple", linewidth=3, linestyle=(0, (5, 10)))
pointdf.plot(ax=ax, markersize=200, color="black")
centerbuffer.boundary.plot(ax=ax, color="black", linewidth=4)
riverdf.plot(ax=ax, color="blue")

plt.show()

enter image description here

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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