无法在图表上显示高度图

Question

无法显示高度图的图形。我在代码中找不到任何问题，但图表是空的。也许有人可以解释如何为这样的项目创建 shp 文件。我不明白：我需要在 shp 中添加一些字段吗？代码有问题？

import matplotlib.pyplot as plt
import geopandas as gpd
import requests
from shapely.geometry import Point
from os.path import join as pth_join
import json
import time

#Configure API request link
request_str = "https://maps.googleapis.com/maps/api/elevation/json?locations=39.7177537,-105.1961914&key=(API here)"

#Define data inputs and outputs
data_root = r"C:\Users\K11El\Desktop\gg"
input_file = "loca.shp"
output_file = input_file[:-4]+"_elev.shp"
#Read input file into a GeoDataFrame
gdf = gpd.read_file(pth_join(data_root,input_file))
gdf_wgs = gdf.to_crs("EPSG:4326")
#Extract x,y from geometry and convert to strings
build_str = gdf_wgs['geometry'].apply(lambda v: ("%s,%s|" % (str(v.y), str(v.x))))
#Limit each request to 512
geom = []
iterations = int(len(build_str)/512) + 1
for x in range(iterations):
    iter_str = build_str[x*512:(x+1)*512]
    #Build concatenated string for API request
    locs = "".join([x for x in iter_str]).rstrip('|')

    #Make request, parse results
    r = requests.get(request_str)
    parsed = json.loads(r.text)

    #Extract elevations from result
    geom = geom + [Point(i['location']['lng'],i['location']
                         ['lat'],i['elevation']) for i in parsed['results']]

    #Slow down API requests so we aren't blocked
    time.sleep(1)
# Save to new file with x,y,z and reproject to CRS of input file
gdf_elev = gpd.GeoDataFrame({'geometry': geom}, crs="EPSG:4326")
gdf_elev.to_crs(gdf.crs).to_file(pth_join(data_root, output_file))
#Check results with a simple plot
gdf_elev['elevation'] = gdf_elev['geometry'].apply(lambda v: (v.z))
gdf_elev.plot("elevation",legend=True,legend_kwds={'label': "Elevation (Meters)"})
plt.show()

Answer 1

使用 osmnx 可以直接使用它从 google

# add x & y columns for osmnx
cities = cities.join(
    cities["geometry"].apply(lambda p: {"x": p.x, "y": p.y}).apply(pd.Series)
)

# create graph
G = ox.utils_graph.graph_from_gdfs(cities, gpd.GeoDataFrame(), graph_attrs=None)
# add elevations using google
G = ox.add_node_elevations_google(G, key)

# now have elevation for all points
gdf = ox.utils_graph.graph_to_gdfs(G, nodes=True, edges=False).set_crs(cities.crs)

样本数据包括高程

osmid	name_left	iso_a3	x	y	elevation	geometry
0	Vatican City	ITA	12.4534	41.9033	24.887	POINT (12.453386544971766 41.903282179960115)
1	San Marino	ITA	12.4418	43.9361	377.216	POINT (12.441770157800141 43.936095834768004)
2	Rome	ITA	12.4813	41.8979	18.527	POINT (12.481312562873995 41.89790148509894)
3	Vaduz	AUT	9.51667	47.1337	513.964	POINT (9.516669472907267 47.13372377429357)
4	Vienna	AUT	16.3647	48.202	185.284	POINT (16.364693096743736 48.20196113681686)
5	Luxembourg	LUX	6.13	49.6117	303.25	POINT (6.130002806227083 49.611660379121076)
6	Monaco	-99	7.40691	43.7396	293.219	POINT (7.406913173465057 43.73964568785249)
7	Andorra	-99	1.51649	42.5	1378.3	POINT (1.5164859605055199 42.5000014435459)
8	Paris	-99	2.33139	48.8686	35.89	POINT (2.33138946713035 48.86863878981461)
9	Ljubljana	SVN	14.515	46.0553	290.773	POINT (14.51496903347413 46.0552883087945)

完整代码

import osmnx as ox
import geopandas as gpd
import pandas as pd
import contextily as ctx

world = gpd.read_file(gpd.datasets.get_path("naturalearth_lowres"))
cities = gpd.read_file(gpd.datasets.get_path("naturalearth_cities")).sjoin(
    world.loc[world["continent"].eq("Europe")]
)

cities = cities.loc[:, ["name_left", "iso_a3", "geometry"]].reset_index(drop=True)
# add x & y columns for osmnx
cities = cities.join(
    cities["geometry"].apply(lambda p: {"x": p.x, "y": p.y}).apply(pd.Series)
)

# create graph
G = ox.utils_graph.graph_from_gdfs(cities, gpd.GeoDataFrame(), graph_attrs=None)
# add elevations using google
G = ox.add_node_elevations_google(G, key)

# now have elevation for all points
gdf = ox.utils_graph.graph_to_gdfs(G, nodes=True, edges=False).set_crs(cities.crs)


ax = gdf.plot(
    column="elevation",
    figsize=[12, 4],
    markersize=100,
    vmax=gdf["elevation"].quantile(0.9),
    legend=True,
)
gdf.apply(
    lambda x: ax.annotate(
        text=x["name_left"],
        xy=x.geometry.coords[0],
        ha="center",
        xytext=(0, 10),
        textcoords="offset points",
    ),
    axis=1,
)

ctx.add_basemap(ax, crs=gdf.crs, source=ctx.providers.Stamen.Terrain)

无法在图表上显示高度图

Cant show heighmap on graph

python

geocoding

geopandas

shapefile

完整代码

输出