【问题标题】:(python) plot 3d surface with colormap as 4th dimension, function of x,y,z(python) 绘制 3d 表面,颜色图为第 4 维,x,y,z 的函数
【发布时间】:2015-12-04 08:30:13
【问题描述】:

我正在尝试绘制一个 3d 表面,其中三个维度中的每一个都位于一个单独的值数组中,并且每个坐标处的表面颜色是 x、y、z 的函数。一种 numpy.pcolormesh 但在 4D 中,而不是 3D。 3D 图由下式给出:

from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm
fig = plt.figure()
ax = fig.gca(projection='3d')
x = np.logspace(-1.,np.log10(5),50)
y = np.linspace(6,9,50)
z = np.linspace(-1,1,50)
colors = LikeBeta(y,range(50),range(50))
ax.plot_trisurf(x,y,z,cmap=colors,linewidth=0.2)

在哪里

def LikeBeta(rho0,r0,beta):
    M0 = 10**rho0*r0_array[r0]**3
    I = cst*M0*sigma_los_beta[beta,:,r0]
    S = dv**2+I
    res = (np.log(S) + (v-u)**2/S).sum()
    return res/2.

cmap=colors 可能是错误的,但问题出在其他地方。我收到以下错误:

----> 8 colors = LikeBeta(y,range(50),range(50))
----> 4     I = cst*M0*sigma_los_beta[beta,:,r0]
    ValueError: operands could not be broadcast together with shapes (50,) (50,353)

确实,sigma_los_beta 是一个我单独评估的数组,其形状为(50,353,50),而这 353 个是我必须拥有的数据。

如何将此函数转换为与plot_trisurf 的其他条目兼容的形式?

抱歉,我无法提供最小的工作代码,因为 dv,v 和 u 是数据。 非常感谢您的帮助。干杯

【问题讨论】:

  • 你可以试试切片吗?您的示例代码中有很多未定义的变量。如果不了解所有内容,就很难提供帮助。

标签: python matplotlib 3d


【解决方案1】:

This 回答解决了 4d 曲面图问题。它使用 matplotlib 的 plot_surface 函数而不是 plot_trisurf

基本上,您希望将 x、y 和 z 变量重塑为相同维度的二维数组。要将第四个维度添加为颜色图,您必须提供另一个与坐标区变量具有相同维度的二维数组。

以下是 3d 绘图的示例代码,其颜色图对应于 x 值。 facecolors 参数用于根据您的喜好更改颜色图。请注意,它的值是从matplotlib.cm.ScalarMappable 类中的to_rgba() 函数获取的。

import matplotlib
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import numpy as np

# domains
x = np.logspace(-1.,np.log10(5),50) # [0.1, 5]
y = np.linspace(6,9,50)             # [6, 9]
z = np.linspace(-1,1,50)            # [-1, 1]

# convert to 2d matrices
Z = np.outer(z.T, z)        # 50x50
X, Y = np.meshgrid(x, y)    # 50x50

# fourth dimention - colormap
# create colormap according to x-value (can use any 50x50 array)
color_dimension = X # change to desired fourth dimension
minn, maxx = color_dimension.min(), color_dimension.max()
norm = matplotlib.colors.Normalize(minn, maxx)
m = plt.cm.ScalarMappable(norm=norm, cmap='jet')
m.set_array([])
fcolors = m.to_rgba(color_dimension)

# plot
fig = plt.figure()
ax = fig.gca(projection='3d')
ax.plot_surface(X,Y,Z, rstride=1, cstride=1, facecolors=fcolors, vmin=minn, vmax=maxx, shade=False)
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')
fig.canvas.show()

我引用的答案(和其他人)提到您应该标准化您的第四维数据。似乎可以通过显式设置颜色图的限制来避免这种情况,就像我在代码示例中所做的那样。

【讨论】:

  • 太棒了!但正如我在问题中所写,我的主要问题在于对第四维数组的评估。按照你的回答和我的代码,我应该写color_dimension = LikeBeta(y,range(50),range(50)),即一个 len(50) 数组。它的元素在LikeBeta 中进行评估,它是 y 的函数以及 x 和 z 的索引,它使用 sigma_los_beta 的组件,一个 (50,353,50) 数组。这导致了我提到的广播问题。一旦我解决了这个问题,我就可以应用你的建议来给表面上色。关于如何解决这个问题的任何想法?谢谢
  • 您能否提供有关LikeBeta 函数中未知数的详细信息?例如。以下类型和维度:r0_array, cst, dv, u, v。如果没有模拟问题的方法,就很难确定问题所在。重要的是该函数应该返回一个与轴变量具有相同维度的数组。另外,您打算执行矩阵乘法还是逐元素乘法?
  • 对不起,你是对的。 r0_array= xcst 只是一个浮点数,dv,u,v 都是 numpy.array 的 len 353(这些是数据点)。因此,为简单起见,您还可以删除这两行代码,只需:def LikeBeta(rho0,r0,beta): M0 = 10**rho0*x[r0]**3 I = M0*sigma_los_beta[beta,:,r0] return I 其中sigma_log_beta 是形状为 (50,353,50) 的 np.array。鉴于我将整个数组 sigma_los_beta[beta,:,r0] 乘以一个标量,它是逐元素的
  • 由于元素乘法中的大小不匹配,您会收到错误消息。逐行跟踪函数中的输出维度,以查看需要在何处应用更改。我不能简单地更改函数来纠正不匹配错误,因为我不知道您想要实现什么。看看numpy.repeat,当你求和时,可能只需要在一个维度上这样做。
  • 我通过使colors 成为具有列表理解的3D 数组,设法解决了数组维度问题。新问题:如何从 3D 数组创建自定义颜色图?显然,上述答案中的 minnmaxx 不适用于 3D 数组。如何从我的 3D 数组创建可映射到颜色图?谢谢
【解决方案2】:

非常感谢 @Frik 的出色 answer,它帮助我实现了 OP 要求的类似情节。

但是,我发现可以对代码进行一些简化,并且可能会引起人们的兴趣。片段和下图。

import matplotlib.pyplot as plt
# This import registers the 3D projection, but is otherwise unused.
from mpl_toolkits.mplot3d import Axes3D  # noqa: F401 unused import
from mpl_toolkits.mplot3d.axes3d import get_test_data
import numpy as np
fig, ax = plt.subplots(subplot_kw={'projection': '3d'})
X, Y, Z = get_test_data(0.05)
C = np.linspace(-5, 5, Z.size).reshape(Z.shape)
scamap = plt.cm.ScalarMappable(cmap='inferno')
fcolors = scamap.to_rgba(C)
ax.plot_surface(X, Y, Z, facecolors=fcolors, cmap='inferno')
fig.colorbar(scamap)
plt.show()

最后,我还想评论一下@Frik 写的内容:

我引用的答案(和其他人)提到您应该标准化您的第四维数据。似乎可以通过像我在代码示例中所做的那样明确设置颜色图的限制来避免这种情况。

我发现这个说法不正确。事实上,如果看一下to_rgba,就会发现有一个norm 关键字,默认设置为True。这正是标准化发生的地方。还包括以下语句:

如果 norm 为 False,则不对输入数据进行归一化,假设在 (0-1) 范围内。

您确实希望您的数据位于 (0-1) 中。

【讨论】:

    【解决方案3】:

    此代码基于 trisurf 演示 http://matplotlib.org/examples/mplot3d/trisurf3d_demo.html

    我添加了一个基于 SO Create own colormap using matplotlib and plot color scale 的函数 make_colormap()

    还添加了一个序列 w=tan(-x*y),它根据该函数生成灰度颜色图。
    您可以使用 cdict 的构造来为其添加更多颜色,但我认为灰度是一个很好的概念证明......

    抱歉,由于缺少最少的工作代码,我无法直接使用您的示例。

    from mpl_toolkits.mplot3d import Axes3D
    from matplotlib import cm
    import matplotlib.pyplot as plt
    import numpy as np
    import matplotlib.colors as mcolors
    
    ###################
    
    def make_colormap(seq):
        """Return a LinearSegmentedColormap
        seq: a sequence of floats and RGB-tuples. The floats should be increasing
        and in the interval (0,1).
        """
        #%
        cdict = {'red': [], 'green': [], 'blue': []}
    
        # make a lin_space with the number of records from seq.     
        x = np.linspace(0,1, len(seq))
        #%
        for i in range(len(seq)):
            segment = x[i]
            tone = seq[i]
            cdict['red'].append([segment, tone, tone])
            cdict['green'].append([segment, tone, tone])
            cdict['blue'].append([segment, tone, tone])
        #%
        return mcolors.LinearSegmentedColormap('CustomMap', cdict)
    
    
    #############################
    
    
    
    n_angles = 36
    n_radii = 8
    
    # An array of radii
    # Does not include radius r=0, this is to eliminate duplicate points
    radii = np.linspace(0.125, 1.0, n_radii)
    
    # An array of angles
    angles = np.linspace(0, 2*np.pi, n_angles, endpoint=False)
    
    # Repeat all angles for each radius
    angles = np.repeat(angles[...,np.newaxis], n_radii, axis=1)
    
    # Convert polar (radii, angles) coords to cartesian (x, y) coords
    # (0, 0) is added here. There are no duplicate points in the (x, y) plane
    x = np.append(0, (radii*np.cos(angles)).flatten())
    y = np.append(0, (radii*np.sin(angles)).flatten())
    
    # Pringle surface
    z = np.sin(-x*y)
    
    
    
    w = np.tan(-x*y)
    colors = make_colormap(w)
    
    
    
    fig = plt.figure()
    ax = fig.gca(projection='3d')
    
    ax.plot_trisurf(x, y, z, cmap=colors, linewidth=0.2)
    
    plt.show()
    

    【讨论】:

    • 非常感谢 PabTorre。但是,如果我尝试在我的代码中评估你的 w 的等价物,这是驱动 3D 表面着色的原因,我会收到广播错误。如果您检查我的LikeBeta 函数,您会注意到我使用三个相同长度的数组调用它。但是,该函数使用不同形状的数组,这会导致广播错误。我无法解决这个问题。你知道我该如何解决这个问题吗?
    • 当我尝试运行你的函数时,r0_array 未定义,该数组的 tge 结构是什么?
    猜你喜欢
    • 2011-04-28
    • 2023-03-17
    • 2015-03-16
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    相关资源
    最近更新 更多