Now that we have been able to create GPX files and a databook for our trip, a logical improvement would be to include water sources. OpenStreetMaps has several tags for identifying water, which seem to have been added over time in an ad-hoc manner. Let us include as many relevant options as we can find.
water_tags = {'natural':['spring','water'],
'water':['pond','lake','reservoir'],
'landuse':['reservoir','basin'],
'amenity':['watering_place','drinking_water','water'],
'waterway':'stream',
'man_made':'water_tank'
}
### ########### FIND WATER SOURCES #############
bbox_gila = (34.35704, 32.77935, -107.65983, -108.94574) # rough box of Gila NF
wgdf = ox.features.features_from_bbox(bbox=bbox_gila, tags=water_tags)
The water features are imported as polygons. That is probably more detail than we need, so let’s reduce the polygon to a point, with the view of being able to import a list of waypoints.
# I really did not want to learn about GeoPandas, but there you go... water_points = wgdf.copy() water_points['geometry'] = water_points['geometry'].centroid
We understand that reducing from polygons to points is an approximation, and might not work so well for rivers and streams and large lakes. Still, we approximate for convenience, and see if the results are close enough to be useful.
We should be able to convert the GeoDataFrame back to a “network” with only nodes and no edges, but our OSMnx library functions do not work well with “networks” with no edges. We could add one dummy edge, or just stay with GeoDataFrames. Eventually we are able to plot results on our map.
Zoom in…
With more than 9k water sources in the Gila imported, this is too many to be practical. Some waypoints seem to be multiple locations along a creek. Other points are so numerous in areas, that we begin to doubt their veracity.
Let us filter out all water sources that are more than X distance from our route.
WJ = WG.copy()
E,D = ox.distance.nearest_edges(JD,X=water_points['x'],Y=water_points['y'],return_dist=True)
max_distance = 2000 #meters
EARTH_RADIUS_M = 6_371_009 # new convenience notation for big numbers
max_distance = max_distance / EARTH_RADIUS_M
# Not clear from documentation, but nearest_edges returns distance in units of
# earth radius if projection is lat/long.
for node,d in zip(WG.nodes(),D):
if d > max_distance:
WJ.remove_node(node)
With the distance filter, our water sources are down to 77. Here they are on the map, and when you hover over a source, it shows some description of the water.
We can also plot the water sources with our route. Water sources appear as darker blue dots.
We could try to add these water sources to our databook, but perhaps the easiest approach is just to import water source waypoints to our favorite GPX app.
Download source code here.
