如何转换contourf值并绘制到3d球体的表面

Question

我正在尝试使用 matplotlibcontourf 生成温度数据的等高线图。然后我想将这些数据映射到 3d 球体上。我正在使用 vpython 来渲染图像。我遇到的问题是多边形不是在球体表面生成的，而且数据中也存在很多间隙。有人可以解释一下如何实现这一点吗？

请注意，我不希望通过从轮廓生成纹理然后用纹理包裹球体来解决这个问题。我希望平移由contourf生成的多边形路径，并将这些路径直接平移到球体上。

from vpython import vector

import numpy as np
from matplotlib import pyplot as plt

from vpython import triangle, vertex, vec


from scipy.interpolate import griddata


def spherical_to_cartesian(lat, lon, radius=3):
    lons = np.radians(lon)
    lats = np.radians(lat)


    x = radius * np.cos(lats) * np.cos(lons)
    y = radius * np.cos(lats) * np.sin(lons)
    z = radius * np.sin(lats)
    return np.array([x, y, z])


shape = (721, 1440)


lats = np.linspace(-90, 90, shape[0])
lons = np.linspace(-180, 180, shape[1])


min_temp = -30
max_temp = 50
temps = np.random.uniform(min_temp, max_temp, size=shape)



lons_grid, lats_grid = np.meshgrid(lons, lats)

new_lons = np.linspace(lons.min(), lons.max(), 72)  # 72 points in longitude
new_lats = np.linspace(lats.min(), lats.max(), 36)  # 36 points in latitude
new_lons_grid, new_lats_grid = np.meshgrid(new_lons, new_lats)

radius = 3

lats_flat = lats_grid.flatten()
lons_flat = lons_grid.flatten()
temps_flat = temps.flatten()

coarse_temps = griddata(
    (lons_flat, lats_flat),
    temps_flat,
    (new_lons_grid, new_lats_grid),
    method='linear'  # Use 'linear', 'nearest', or 'cubic' as needed
)


norm = plt.Normalize(coarse_temps.min(), vmax=coarse_temps.max())

cmap = plt.get_cmap("inferno", 100)

levels = 100

levels = np.linspace(coarse_temps.min(), coarse_temps.max(), levels + 1)

bucketed_data = np.digitize(coarse_temps, levels) - 1

fig, ax = plt.subplots()
contours = ax.contourf(new_lons_grid, new_lats_grid, bucketed_data, levels=np.arange(len(levels)))
plt.close(fig)



def create_polygon(region, temp_color):
    if len(region) < 3:
        # Can't form a polygon with fewer than 3 points
        return

    # Calculate the centroid of the region
    centroid = vec(sum(p[0] for p in region) / len(region),
                   sum(p[1] for p in region) / len(region),
                   sum(p[2] for p in region) / len(region))

    # Tessellate the region into triangles
    for i in range(len(region) - 1):
        v0 = vec(*region[i])
        v1 = vec(*region[i + 1])
        triangle(
            v0=vertex(pos=v0, color=temp_color),
            v1=vertex(pos=v1, color=temp_color),
            v2=vertex(pos=centroid, color=temp_color)
        )


    v0 = vec(*region[-1])
    v1 = vec(*region[0])
    triangle(
        v0=vertex(pos=v0, color=temp_color),
        v1=vertex(pos=v1, color=temp_color),
        v2=vertex(pos=centroid, color=temp_color)
    )


def split_contours(segs, kinds=None):
    if kinds is None:
        return segs  # nothing to be done
    new_segs = []
    for i, seg in enumerate(segs):
        segkinds = kinds[i]
        boundaries = [0] + list(np.nonzero(segkinds == 79)[0])
        for b in range(len(boundaries) - 1):
            new_segs.append(seg[boundaries[b] + (1 if b > 0 else 0):boundaries[b + 1]])
    return new_segs



allsegs = contours.allsegs
allkinds = contours.allkinds
colors = cmap(norm(coarse_temps))
rgb_colors = [
    tuple(int(c * 255) for c in color[:3])
    for color in colors.reshape(-1, colors.shape[-1])
]
for clev in range(len(contours.allsegs)):
    kinds = None if allkinds is None else allkinds[clev]
    segs = split_contours(allsegs[clev], kinds)
    rgb = cmap(clev)[:3]
    coords = [spherical_to_cartesian(lat, lon, radius=3) for seg in segs for lon, lat in seg]
    temp_color = vector(*rgb)
    create_polygon(coords, temp_color)



import time
while True:
    time.sleep(0.03)

这就是球体的渲染方式。我希望将轮廓渲染在球体的表面上，就像 2d 轮廓投影到 3d 表面上一样。

Answer 1

感谢 Reddit 上的

misho88

，他发布了 this 并建议查看使用

matplotlib.tri.Triangulation

和

tricontourf

以下作品：

import numpy as np
import matplotlib.tri as tri
from vpython import vector, triangle, vertex, vec
import matplotlib.pyplot as plt


# Function to convert spherical to Cartesian coordinates
def spherical_to_cartesian(lat, lon, radius=3):
    lons = np.radians(lon)
    lats = np.radians(lat)
    x = radius * np.cos(lats) * np.cos(lons)
    y = radius * np.cos(lats) * np.sin(lons)
    z = radius * np.sin(lats)
    return np.array([x, y, z])

shape = (721, 1440)


lats = np.linspace(-90, 90, shape[0])
lons = np.linspace(-180, 180, shape[1])


min_temp = -30
max_temp = 50
temps = np.random.uniform(min_temp, max_temp, size=shape)

lons_grid, lats_grid = np.meshgrid(lons, lats)

new_lons = np.linspace(lons.min(), lons.max(), 72)
new_lats = np.linspace(lats.min(), lats.max(), 36)
new_lons_grid, new_lats_grid = np.meshgrid(new_lons, new_lats)

radius = 3

lats_flat = lats_grid.flatten()
lons_flat = lons_grid.flatten()
temps_flat = temps.flatten()

# Randomly sample a subset of the points
num_samples = 10000  # Adjust this number to control triangle density
indices = np.random.choice(len(lats_flat), size=num_samples, replace=False)
lats_sampled = lats_flat[indices]
lons_sampled = lons_flat[indices]
temps_sampled = temps_flat[indices]

triang = tri.Triangulation(lons_sampled, lats_sampled)


fig, ax = plt.subplots()
contour = ax.tricontourf(triang, temps_sampled, levels=100, cmap='inferno')
plt.colorbar(contour)
plt.close(fig)

for tri_indices in triang.triangles:

    vertices = []
    for idx in tri_indices:
        xi, yi, zi = spherical_to_cartesian(lats_sampled[idx], lons_sampled[idx])
        temp = temps_sampled[idx]
        color_value = plt.cm.inferno((temp - min_temp) / (max_temp - min_temp))
        vertices.append(vertex(pos=vector(xi, yi, zi), color=vec(*color_value[:3])))
    triangle(v0=vertices[0], v1=vertices[1], v2=vertices[2])


import time
while True:
    time.sleep(0.03)