如何沿给定线从栅格中提取价值概况？

Question

如何沿 Python 中给定的 shapefile 线从栅格中提取值的配置文件？

我正在努力寻找一种方法来从栅格 (geotiff) 中提取值的轮廓（例如地形轮廓）。来自基于多边形的栅格的库 Rasterio has a method to clip/extract 值，但我找不到线形文件的等效方法。

有一个 basic method with scipy，但它并不像 rasterio 可以提供的基于更高级别工具箱的方法那样固有地保存地理信息。

换句话说，我在 Python 中寻找与 QGIS 中的 Terrain Profile 工具提供的等效项。

谢谢

Answer 1

这与提取多边形有点不同，因为您希望按照触摸的顺序对线条触摸的每个像素进行采样（多边形方法不关心像素顺序）。

看起来可以调整 this approach to use rasterio instead. Given a line read from a shapefile using geopandas or fiona as a shapely object, you use the endpoints to derive a new equidistant projection that you use as dst_crs in a WarpedVRT 并从中读取像素值。看起来您需要根据要采样的像素数来计算线条的长度，这是 WarpedVRT.

的宽度参数

如果您的线不是端点之间的近似直线，则可能需要进一步调整此方法。

如果您只想获取线下的原始像素值，您应该能够为每条线直接使用 mask in rasterio or rasterize。在行的情况下，您可能希望使用 all_touched=True。

Answer 2

我遇到了类似的问题并找到了适合我的解决方案。该解决方案使用 shapely 对 line/lines 上的点进行采样，然后从 GeoTiff 访问相应的值，因此提取的配置文件遵循线的方向。这是我最终得到的方法：

def extract_along_line(xarr, line, n_samples=256):
    profile = []

    for i in range(n_samples):
        # get next point on the line
        point = line.interpolate(i / n_samples - 1., normalized=True)
        # access the nearest pixel in the xarray
        value = xarr.sel(x=point.x, y=point.y, method="nearest").data
        profile.append(value)
        
    return profile

这是一个工作示例，数据来自 copernicus-dem 数据库，直线是接收到的图块的对角线：

import rioxarray
import shapely.geometry
import matplotlib.pyplot as plt

sample_tif = ('https://elevationeuwest.blob.core.windows.net/copernicus-dem/'
              'COP30_hh/Copernicus_DSM_COG_10_N35_00_E138_00_DEM.tif')

# Load xarray
tile = rioxarray.open_rasterio(sample_tif).squeeze()
# create a line (here its the diagonal of tile)
line = shapely.geometry.MultiLineString([[
            [tile.x[-1],tile.y[-1]],
            [tile.x[0], tile.y[0]]]])

# use the method from above to extract the profile
profile = extract_along_line(tile, line)
plt.plot(profile)
plt.show()

如何沿给定线从栅格中提取价值概况？

How to extract a profile of value from a raster along a given line?

python

fiona

geopandas

rasterio