如何在 Python 中将等值线图与阴影栅格相结合？

Question

我想在地图上绘制区域特征，但人口密度非常不均匀，较大的图块会误导人们注意。想想邮政编码的平均值（比如考试成绩）。

高分辨率地图可用于区分有人居住的地点，甚至是其中的密度。下面的 Python 代码确实会根据每个像素的平均此类密度生成彩色光栅。

然而，我真正需要的是从同一区域（在本例中为匈牙利的邮政编码）的等值线图着色，但着色只影响无论如何都会出现在光栅上的点。光栅只能确定像素的伽马（或者可能是某些 3D 模拟中的高度）。解决这个问题的好方法是什么？

一个rasterio.mask.mask不知何故？

（顺便说一句，用邮政编码边界覆盖也很好，但我更好地理解了如何使用 GeoViews。）

import rasterio
import os
import datashader as ds
from datashader import transfer_functions as tf
import xarray as xr
from matplotlib.cm import viridis

# download a GeoTIFF from this location: https://data.humdata.org/dataset/hungary-high-resolution-population-density-maps-demographic-estimates
data_path = '~/Downloads/'
file_name = 'HUN_youth_15_24.tif'  # young people
file_path = os.path.join(data_path, file_name)
src = rasterio.open(file_path)
da = xr.open_rasterio(file_path)
cvs = ds.Canvas(plot_width=5120, plot_height=2880)
img = tf.shade(cvs.raster(da,layer=1), cmap=viridis)
ds.utils.export_image(img, "map", export_path=data_path, fmt=".png")

Answer 1

Datashader 可让您将多种类型的数据组合成一个常见的栅格形状，您可以在其中使用基于 NumPy 的 xarray 操作进行任何制作或过滤。例如。您可以将等值线渲染为多边形，然后屏蔽掉无人居住的区域。如何按区域归一化取决于您，并且可能会变得非常复杂，但是一旦您准确定义了您打算做什么就应该可行。有关如何执行此操作的示例，请参阅 https://examples.pyviz.org/nyc_taxi/nyc_taxi.html 处的 transform 代码，如：

def transform(overlay):
    picks = overlay.get(0).redim(pickup_x='x', pickup_y='y')
    drops = overlay.get(1).redim(dropoff_x='x', dropoff_y='y')
    pick_agg = picks.data.Count.data
    drop_agg = drops.data.Count.data
    more_picks = picks.clone(picks.data.where(pick_agg>drop_agg))
    more_drops = drops.clone(drops.data.where(drop_agg>pick_agg))
    return (hd.shade(more_drops, cmap=['lightcyan', "blue"]) *
            hd.shade(more_picks, cmap=['mistyrose', "red"]))

picks = hv.Points(df, ['pickup_x',  'pickup_y'])
drops = hv.Points(df, ['dropoff_x', 'dropoff_y'])
((hd.rasterize(picks) * hd.rasterize(drops))).apply(transform).opts(
    bgcolor='white', xaxis=None, yaxis=None, width=900, height=500)

这里并没有真正屏蔽任何东西，但希望你能看到屏蔽是如何工作的；只需获取一些栅格化对象，然后使用其他一些栅格化对象进行数学运算。这里的所有步骤都是在一个使用 HoloViews 对象的函数中完成的，这样您就可以拥有一个实时的交互式绘图，但是您可能希望使用 datashader.org 处的更基本的代码来解决该方法，您只需要处理xarray 对象而不是 HoloViews 管道；然后，您可以将对单个 xarray 所做的操作转换为 HoloViews 管道，然后允许与平移、缩放、轴等进行完全交互使用。

Answer 2

我不确定我是否理解，所以如果我弄错了请告诉我。如果我理解得很好，你可以只使用 numpy 来实现你想要的（我相信将它翻译成 xarray 会很容易）：

# ---- snipped code already in the question -----
import numpy as np
import matplotlib.pyplot as plt

#  fake a choropleth in a dirty, fast way
height, width = 2880, 5120
choropleth = np.empty((height, width, 3,), dtype=np.uint8)
CHUNKS = 10
x_size = width // CHUNKS
for x_step, x in enumerate(range(0, width, width // CHUNKS)):
    y_size = height // CHUNKS
    for y_step, y in enumerate(range(0, height, height // CHUNKS)):
        choropleth[y: y+y_size, x: x+x_size] = (255-x_step*255//CHUNKS,
                                                0, y_step*255//CHUNKS)
plt.figure("Fake Choropleth")
plt.imshow(choropleth)

# Option 1: play with alpha only
outimage = np.empty((height, width, 4,), dtype=np.uint8)  # RGBA image
outimage[:, :, 3] = img  # Set alpha channel
outimage[:, :, :3] = choropleth  # Set color
plt.figure("Alpha filter only")
plt.imshow(outimage)

# Option 2: clear the empty points
outimage[img == 0, :3] = 0  # White. use 0 for black
plt.figure("Points erased")
plt.imshow(outimage[:,:,:3])  # change to 'outimage' to see the image with alpha

结果：虚拟 choroplet

Alpha 过滤图

黑色背景，无 alpha 滤镜请注意，由于 matplotlib 的抗锯齿功能，图像可能看起来有所不同。

如何在 Python 中将等值线图与阴影栅格相结合？

How can a choropleth map be combined with a shaded raster in Python?

python

rasterio

datashader

geoviews