You need to enable JavaScript to run this app.
优惠活动
大模型
产品
解决方案
定价
更多

使用Cartopy添加几何图形时,对数归一化值显示颜色过浅

斯堪的纳维亚市政区火灾强度可视化异常问题

我尝试用contourf在斯堪的纳维亚地图上按市政区(.shp文件)绘制.csv格式的火灾强度数据。由于数据值范围跨度极大(最大值达10³量级),我采用对数尺度展示低值区间细节。之前代码运行正常,但现在执行相同代码时,地图颜色过浅,既不匹配实际数据也不符合色条显示,推测和库版本更新有关,但找不到修复方法。另外,也想了解更简便的对数尺度展示方法。

我的代码

import cartopy.crs as ccrs
import numpy as np
from shapely.geometry import Polygon
import cartopy
from cartopy.io import shapereader
import matplotlib
import pandas
import geopandas
import cartopy.io.img_tiles as cimgt
import matplotlib.pyplot as plt
import cartopy.feature as feature
import shapely

projection = ccrs.Mercator()
transform = ccrs.PlateCarree()
extent = [4.5, 31.9, 54.7, 71.24]
cmap = matplotlib.cm.get_cmap('Reds')

# Get borders
resolution = '10m'
category = 'cultural'
name = 'admin_0_countries'

shpfilename = shapereader.natural_earth(resolution, category, name)
shpfilenamemunicipalities = '/path/to/ScandinaviaKommuner.shp'

# Read the shapefile using geopandas
df = geopandas.read_file(shpfilename)
df2 = geopandas.read_file(shpfilenamemunicipalities)

# Select the countries of interest
scan3 = df[ df['ADMIN'].isin(['Norway', 'Finland', 'Sweden']) ]
scan3_dissolved = scan3.dissolve(by='LEVEL')
poly = [scan3_dissolved['geometry'].values[0]]

# Create the mask
stamen_terrain = cimgt.Stamen('terrain_background')
st_proj = stamen_terrain.crs  #projection used by Stamen images
ll_proj = ccrs.PlateCarree()

def rect_from_bound(xmin, xmax, ymin, ymax):
    """Returns list of (x,y)'s for a rectangle"""
    xs = [xmax, xmin, xmin, xmax, xmax]
    ys = [ymax, ymax, ymin, ymin, ymax]
    return [(x, y) for x, y in zip(xs, ys)]

pad1 = .1  #padding, degrees unit
exts = [poly[0].bounds[0] - pad1, poly[0].bounds[2] + pad1, poly[0].bounds[1] - pad1, poly[0].bounds[3] + pad1];
msk = Polygon(rect_from_bound(*exts)).difference(poly[0].simplify(0.01))
msk_stm  = st_proj.project_geometry(msk, ll_proj)  # project geometry to the projection used by stamen

# Load the .csv data
firedata = pandas.read_csv('/path/to/countrylevel_firedata.csv', delimiter=';', decimal='.')
firedatamunicipalities = pandas.read_csv('/path/to/municipalitylevel_firedata.csv', delimiter=';', decimal='.')

# Subset .csv data
countries = firedata['Country'].tolist()
municipalities = firedatamunicipalities['Municipality'].tolist()
fireextents = firedatamunicipalities['FireIntensity'].tolist()

# Normalise the fire data to between 0 and 1 to extract the colour and log transform
fireextents_log = list(range(0,966))
for fireextent, log in zip(fireextents, fireextents_log):
    if fireextent>0:
        fireextents_log[log] = np.log(fireextent)
    else:
        fireextents_log[log] = np.nan

fireextents_norm = (fireextents_log-np.nanmin(fireextents_log))/(np.nanmax(fireextents_log) - np.nanmin(fireextents_log))

# Create a plot
fig, axs = plt.subplots(nrows = 1, ncols = 1,
                        subplot_kw = {'projection': projection},
                        layout='constrained')

axs.set_extent(extent)
axs.add_feature(feature.BORDERS, linewidth=2.5, zorder=15)
axs.set_extent(extent)
axs.axis('off')
axs.add_geometries(shapely.get_parts(poly[0]), crs=transform, facecolor='none',
                      edgecolor='black', linewidth=2.5, zorder=18)
axs.add_geometries(msk_stm, st_proj, facecolor='#f3f8ff', edgecolor='none',
                      zorder=16)

for municipality, fireextent_norm in zip(municipalities, fireextents_norm):
    poly3 = df2.loc[df2['Municipali'] == municipality]['geometry'].values
    if fireextent_norm > 0:
        rgba2 = cmap(fireextent_norm)
    else:
        rgba2 = 'white'
    axs.add_geometries(poly3, crs=transform, facecolor=rgba2, edgecolor='black', linewidth=1, zorder=14)

dummy_scat = axs.scatter(fireextents, fireextents, c=fireextents, cmap=cmap, zorder=0,
                          norm=matplotlib.colors.LogNorm())
cbar = fig.colorbar(dummy_scat, ax=axs, orientation='horizontal', shrink=1, pad=0.01)
cbar.set_label('Fire intensity (ha burnt forest/fire) \nover the period 2008-2024', size=65, labelpad=15)
cbar.ax.tick_params(labelsize=60)

plt.show()

部分火灾强度数据示例

print(fireextents)
>>> [4183.6667, 1378.0, 1065.0, 991.6, 908.0, 612.0, 548.0625, 520.0, 412.7778, 331.0, 324.0, 280.5, 259.3333, 259.0, 234.0, 222.0, 212.6667, 199.0, 169.0, 158.3333, 153.5, 148.6, 147.5, 145.0, 133.8, 130.0, 108.3333, 105.25, 103.0, 103.0, 90.0, 86.0, 84.3, 82.1667, 80.6, 80.0, 79.8571, 75.0, 74.0, 71.0, 65.3, 63.0, 59.6, 58.0, 58.0, 57.6667, 55.0, 55.0, 54.9, 54.8636, 53.25, 51.0, 50.4667, 49.0, 49.0, 49.0, 48.7143, 48.0, 47.1071, 47.0, 46.5, 46.0, 45.0, 43.0, 42.2, 42.0, 41.0, 40.0, 38.2, 35.4, 34.9, 34.9, 34.5, 34.0, 34.0, 34.0, 33.8, 33.0, 33.0, 33.0, 32.8, 32.0, 31.0, 30.3333, 30.0, 30.0, 30.0, 29.6667, 29.5, 29.4, 29.3, 29.2857, 28.1667, 28.0, 28.0, 28.0, 27.0, 27.0, 26.0, 26.0, 26.0, 26.0, 26.0, 26.0, 25.8, 25.0, 25.0, 25.0, 25.0, 24.6, 24.5, 24.0, 24.0, 23.5, 23.1429, 23.0, 23.0, 22.0, 21.3, 19.0, 19.0, 19.0, 19.0, 19.0, 18.0, 17.5, 17.0, 16.5, 15.9167, 15.8, 15.5, 15.5, 15.3, 15.0714, 15.0, 15.0, 15.0, 14.5, 14.0, 13.5, 13.3333, 13.0, 12.5, 12.5, 12.2941, 12.0588, 12.0, 12.0, 12.0, 11.6, 11.0, 10.6667, 10.0, 9.6667, 9.6667, 9.5, 9.5, 9.0, 8.5714, 8.5, 8.3, 7.1667, 7.0, 7.0, 7.0, 6.6, 6.0, 6.0, 6.0, 6.0, 5.5, 5.0, 5.0, 5.0, 5.0, 5.0, 5.0, 5.0, 5.0, 5.0, 4.5, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 3.75, 3.5, 3.5, 3.4, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 2.5, 2.2, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 1.5, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, ... 0.0, 0.0, 0.0]

解决方案

1. 颜色异常的核心原因

你手动做的对数归一化和色条使用的LogNorm逻辑完全脱节:

  • 代码中先对数据取对数再线性归一化到[0,1],用这个值映射颜色
  • 但色条直接用LogNorm对原始数据做对数映射,两者的颜色映射规则不一致,库版本更新后matplotlib的colormap细节变化,导致颜色显示异常。

2. 修复代码

统一使用LogNorm处理颜色映射,去掉手动的对数归一化步骤:

替换数据处理代码

# 删除原有的fireextents_log和fireextents_norm代码,替换为:
# 过滤掉0值,计算LogNorm的范围
valid_fire = [x for x in fireextents if x > 0]
norm = matplotlib.colors.LogNorm(vmin=np.min(valid_fire), vmax=np.max(valid_fire))

修改市政区绘制循环

直接用norm对原始火灾值做映射,确保和色条逻辑一致:

for municipality, fire_val in zip(municipalities, fireextents):
    poly3 = df2.loc[df2['Municipali'] == municipality]['geometry'].values
    if fire_val > 0:
        rgba2 = cmap(norm(fire_val))
    else:
        rgba2 = 'white'
    axs.add_geometries(poly3, crs=transform, facecolor=rgba2, edgecolor='black', linewidth=1, zorder=14)

简化色条创建

不需要dummy散点,直接用ScalarMappable创建和映射逻辑一致的色条:

# 删除dummy_scat相关代码,替换为:
sm = matplotlib.cm.ScalarMappable(norm=norm, cmap=cmap)
sm.set_array([])  # 仅用于创建色条,不需要数据
cbar = fig.colorbar(sm, ax=axs, orientation='horizontal', shrink=1, pad=0.01)
cbar.set_label('Fire intensity (ha burnt forest/fire) \nover the period 2008-2024', size=65, labelpad=15)
cbar.ax.tick_params(labelsize=60)

3. 更简便的对数尺度展示方法

用GeoPandas的合并+绘图功能,一步完成,避免手动循环的繁琐和错误:

# 合并市政区矢量数据和火灾数据
merged_df = df2.merge(firedatamunicipalities, left_on='Municipali', right_on='Municipality')

# 直接绘制,自动处理颜色映射和0值(缺失值)
merged_df.plot(
    column='FireIntensity', 
    ax=axs, 
    transform=transform, 
    cmap='Reds', 
    norm=matplotlib.colors.LogNorm(),
    edgecolor='black', 
    linewidth=1, 
    zorder=14,
    missing_kwds={'color': 'white'}  # 0值显示为白色
)

这种方法不仅代码更简洁,还能自动处理数据匹配,避免手动循环可能出现的索引错误。

额外优化:支持0值的对数尺度

如果需要展示0值的分布,而不是设为白色,可以用SymLogNorm,它支持包含0的对数尺度:

norm = matplotlib.colors.SymLogNorm(linthresh=1, vmin=0, vmax=np.max(valid_fire), base=10)

linthresh参数指定线性映射的阈值,小于该值的部分用线性映射,大于的部分用对数映射,既保留对数尺度的细节,又能展示0值。

内容的提问来源于stack exchange,提问作者Elsri

相关产品推荐
方舟 Agent Plan

超全模态模型 × Harness 升级,最新支持 Deepseek-V4.1-Flash、GLM-5.3 系列、Doubao-Seedream-5.0-pro、Kimi-K3 (部分), 限时 9.9 元起

最近更新时间:2026.06.18 04:29:51