在 Matplotlib 中根据像素值设置透明度

Question

我正在尝试使用 matplotlib 为我正在写的论文绘制一些图形。我在 2D numpy 数组中有两组数据：一个 ascii 山体阴影栅格，我可以愉快地使用它来绘制和调整：

import matplotlib.pyplot as pp
import numpy as np

hillshade = np.genfromtxt('hs.asc', delimiter=' ', skip_header=6)[:,:-1]

pp.imshow(hillshade, vmin=0, vmax=255)
pp.gray()
pp.show()

这给出了：

第二个 ascii 栅格描绘了流经景观的河流的属性。该数据可以按照与上述相同的方式绘制，但是数组中与河流网络不对应的值被指定为无数据值 -9999。目的是将无数据值设置为透明，以便河流值覆盖山体阴影。

这是河流数据，理想情况下这里表示为 0 的每个像素都是完全透明的。

对此进行了一些研究后，我似乎可以将数据转换为 RGBA 数组并设置 alpha 值以仅使不需要的单元格透明。但是，河流数组中的值是浮点数，无法转换（因为原始值是图形的整个点），并且我相信如果使用 RGBA 格式，

imshow

函数只能采用无符号整数。

有什么办法可以绕过这个限制吗？我曾希望我可以简单地创建一个包含像素值和 alpha 值的元组并像这样绘制它们，但这似乎不可能。

我还尝试过使用

PIL

来尝试创建一个没有数据值透明的河流数据的PNG文件，但这似乎会自动将像素值缩放到0-255，从而丢失我需要的值保存。

我欢迎任何人对这个问题有任何见解。

Answer 1

只需遮罩你的“河流”数组。

例如

rivers = np.ma.masked_where(rivers == 0, rivers)

作为以这种方式叠加两个图的快速示例：

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.cm as cm

# Generate some data...
gray_data = np.arange(10000).reshape(100, 100)

masked_data = np.random.random((100,100))
masked_data = np.ma.masked_where(masked_data < 0.9, masked_data)

# Overlay the two images
fig, ax = plt.subplots()
ax.imshow(gray_data, cmap=cm.gray)
ax.imshow(masked_data, cmap=cm.jet, interpolation='none')
plt.show()

另外，顺便说一句，

imshow

很乐意接受 RGBA 格式的浮点数。它只是期望一切都在 0 到 1 之间。

Answer 2

在不使用屏蔽数组的情况下执行此操作的另一种方法是设置颜色图如何处理低于

clim

最小值的剪切值（无耻地使用 Joe Kington 的示例）：

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.cm as cm

# Generate some data...
gray_data = np.arange(10000).reshape(100, 100)

masked_data = np.random.random((100,100))

my_cmap = cm.jet
my_cmap.set_under('k', alpha=0)


# Overlay the two images
fig, ax = plt.subplots()
ax.imshow(gray_data, cmap=cm.gray)
im = ax.imshow(masked_data, cmap=my_cmap, 
          interpolation='none', 
          clim=[0.9, 1])
plt.show()

还有一个

set_over

用于剪掉顶部，还有一个

set_bad

用于设置颜色图如何处理数据中的“坏”值。

这样做的一个优点是您可以通过调整

clim

和

im.set_clim([bot, top])

来更改阈值

Answer 3

另一个选项是设置所有对

np.nan

保持透明的单元格（不确定这里什么更有效，我猜tacaswell的答案基于

clim

将是最快的）。改编示例Joe Kington 的答案：

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.cm as cm

# Generate some data...
gray_data = np.arange(10000).reshape(100, 100)

masked_data = np.random.random((100,100))
masked_data[np.where(masked_data < 0.9)] = np.nan

# Overlay the two images
fig, ax = plt.subplots()
ax.imshow(gray_data, cmap=cm.gray)
ax.imshow(masked_data, cmap=cm.jet, interpolation='none')
plt.show()

请注意，对于

dtype=bool

数组，您不应为了 PEP 8 (E712) 而遵循 IDE 的建议来比较

masked_data is True

，而应坚持使用

masked_data == True

进行逐元素比较，否则屏蔽将失败：

Answer 4

另一个策略是输出河流支流网络。这通常是水文流量演算工具中的一个选项。我不确定您在这里使用的是哪个工具，但下面是使用 Python

pysheds

库的示例，它需要栅格类型文件，因此使用 tiff 文件而不是 2D numpy 数组：

from pysheds.grid import Grid
import matplotlib.pyplot as plt
import numpy as np

# Instantiate grid from raster
grid = Grid.from_raster('test1.tif')
dem = grid.read_raster('test1.tif')

# Fill pits
pit_filled_dem = grid.fill_pits(dem)

# Fill depressions
flooded_dem = grid.fill_depressions(pit_filled_dem)

# Resolve flats and compute flow directions
inflated_dem = grid.resolve_flats(flooded_dem)
fdir = grid.flowdir(inflated_dem)

# Compute accumulation
acc = grid.accumulation(fdir)

# Extract river branch network for accumulations > 1000 units
branches = grid.extract_river_network(fdir, acc > 1000)

# Create fig and axes objects of matplotlib figure
fig, ax = plt.subplots()

# Set limits and aspect and tick formatting
plt.xlim(grid.bbox[0], grid.bbox[2])
plt.ylim(grid.bbox[1], grid.bbox[3])
ax.set_aspect('equal')
ax.ticklabel_format(style='sci', scilimits=(0,0), axis='both')

# Set the axes color to gray to demonstrate non-river pixels 
ax.set_facecolor('lightgray')

# Plot the river branch network using for loop
for branch in branches['features']:
    line = np.asarray(branch['geometry']['coordinates'])
    plt.plot(line[:, 0], line[:, 1], color='k')

# You can also export the river branch network to a geodataframe
for branch in branches['features']:
    line_coords = branch['geometry']['coordinates']
    line = LineString(line_coords)
    line_strings.append(line)

gdf = gpd.GeoDataFrame(geometry=line_strings)
gdf.plot(ax=ax)