备注
互动在线版： 
使用GeoPandas和Rasterio采样点数据#
此示例说明如何将GeoPandas与Rasterio一起使用。 Rasterio 是用于读取和写入栅格数据的包。
在本例中，使用一组矢量点对这些点处的栅格数据进行采样。
所使用的栅格数据是用于哨兵数据的哥白尼哨兵数据2018。
[1]:
import geopandas
import rasterio
import matplotlib.pyplot as plt
from shapely.geometry import Point
---------------------------------------------------------------------------
ModuleNotFoundError Traceback (most recent call last)
Input In [1], in <cell line: 2>()
1 import geopandas
----> 2 import rasterio
3 import matplotlib.pyplot as plt
4 from shapely.geometry import Point
ModuleNotFoundError: No module named 'rasterio'
创建示例矢量数据#
根据一组点生成地理数据框
[2]:
# Create sampling points
points = [Point(625466, 5621289), Point(626082, 5621627), Point(627116, 5621680), Point(625095, 5622358)]
gdf = geopandas.GeoDataFrame([1, 2, 3, 4], geometry=points, crs=32630)
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
Input In [2], in <cell line: 2>()
1 # Create sampling points
----> 2 points = [Point(625466, 5621289), Point(626082, 5621627), Point(627116, 5621680), Point(625095, 5622358)]
3 gdf = geopandas.GeoDataFrame([1, 2, 3, 4], geometry=points, crs=32630)
NameError: name 'Point' is not defined
这个 GeoDataFrame 如下所示：
[3]:
gdf.head()
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
Input In [3], in <cell line: 1>()
----> 1 gdf.head()
NameError: name 'gdf' is not defined
打开栅格数据#
使用 rasterio 打开要采样的栅格数据的步骤
[4]:
src = rasterio.open('s2a_l2a_fishbourne.tif')
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
Input In [4], in <cell line: 1>()
----> 1 src = rasterio.open('s2a_l2a_fishbourne.tif')
NameError: name 'rasterio' is not defined
让我们来看看叠加点数据的栅格数据。
[5]:
from rasterio.plot import show
fig, ax = plt.subplots()
# transform rasterio plot to real world coords
extent=[src.bounds[0], src.bounds[2], src.bounds[1], src.bounds[3]]
ax = rasterio.plot.show(src, extent=extent, ax=ax, cmap='pink')
gdf.plot(ax=ax)
---------------------------------------------------------------------------
ModuleNotFoundError Traceback (most recent call last)
Input In [5], in <cell line: 1>()
----> 1 from rasterio.plot import show
3 fig, ax = plt.subplots()
5 # transform rasterio plot to real world coords
ModuleNotFoundError: No module named 'rasterio'
采样数据#
Rasterio需要x，y格式的坐标列表，而不是geomentry列中的点。
这可以使用下面的代码来实现
[6]:
coord_list = [(x,y) for x,y in zip(gdf['geometry'].x , gdf['geometry'].y)]
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
Input In [6], in <cell line: 1>()
----> 1 coord_list = [(x,y) for x,y in zip(gdf['geometry'].x , gdf['geometry'].y)]
NameError: name 'gdf' is not defined
对数据进行采样，并将结果存储在一个名为 value 。请注意，如果图像有多个波段，则会为每个波段返回值。
[7]:
gdf['value'] = [x for x in src.sample(coord_list)]
gdf.head()
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
Input In [7], in <cell line: 1>()
----> 1 gdf['value'] = [x for x in src.sample(coord_list)]
2 gdf.head()
NameError: name 'src' is not defined