如何将线向量转换为线密度栅格？

Question

我正在尝试将多条线的矢量图层转换为像素值等于穿过该单元格的线数的栅格。

我不知道从哪里开始。任何帮助，将不胜感激。谢谢！

Answer 1

import geopandas as gpd
import numpy as np
import pandas as pd
from shapely.geometry import Polygon
from shapely.geometry import LineString
from geopandas import sjoin
import rasterio
import xarray as xr
from geocube.api.core import make_geocube
import os

# make some linestrings
line1 = LineString([(17.45, 51.11), (17.46, 51.35)])
line2 = LineString([(17.2, 51.8), (17.23, 51.72)])
line3 = LineString([(16.91, 51.3), (16.84, 50.65)])
line4 = LineString([(16.31, 51.44), (16.37, 51.4)])
line5 = LineString([(16.38, 51.42), (16.43, 51.49)])

# make linestring gdf and save to .shp
line_gdf = gpd.GeoDataFrame([line1, line2, line3, line4, line5],
                            geometry=[line1, line2, line3, line4, line5])

line_gdf.crs = "EPSG:4326"
line_gdf = line_gdf[['geometry']]
line_gdf.to_file('line_gdf.shp')

# make shapely polygon grid 
xmin,ymin,xmax,ymax = [16, 50, 18, 52]

length = 0.25
width = 0.25

cols = list(np.arange(np.floor(xmin), np.ceil(xmax), width))
rows = list(np.arange(np.floor(ymin), np.ceil(ymax), length))
rows.reverse()

polygons = []
for x in cols:
    for y in rows:
        polygons.append( Polygon([(x,y), (x+width, y), (x+width, y-length), (x, y-length)]) )

#make grid_gdf, save to .shp
grid_gdf = gpd.GeoDataFrame({'geometry':polygons})
grid_gdf.crs = "EPSG:4326"
grid_gdf.to_file('grid_gdf.shp')

# loop through grid cells and count how many lines intersect each grid cell
int_count_dict = {}
for idx, row in grid_gdf.iterrows():
    lines_in_polygons_gdf = gpd.sjoin(line_gdf,
                                      grid_gdf.iloc[idx:idx+1,:],
                                      op='intersects',
                                      how='inner')
    # store line intersections as dict values
    int_count_dict[idx] = lines_in_polygons_gdf.shape[0]

#map intersections to grid_gdf into count column
grid_gdf['count'] = grid_gdf.index.map(int_count_dict)

#make xarray of grid_gdf
cube = make_geocube(vector_data=grid_gdf, measurements=['count'], resolution=(0.25, -0.25))

#export to .tif
cube["count"].rio.to_raster("line_count_raster.tif")

print('gpd version:', gpd.__version__)
print('np version:', np.__version__)
print('pd version:', pd.__version__)
print('shapely version:', shapely.__version__)
print('rasterio version:', rasterio.__version__)
print('xarray version:', xr.__version__)
print('geocube version:', geocube.__version__)

gpd version: 0.7.0
np version: 1.18.1
pd version: 1.0.3
shapely version: 1.7.0
rasterio version: 1.1.3
xarray version: 0.15.1
geocube version: 0.0.11

如何将线向量转换为线密度栅格？

How can line vectors be converted into a raster of line densities?

python

gis

spatial

geopandas

rasterio