将图像映射到球体上并绘制 3D 轨迹
Map an image onto a sphere and plot 3D trajectories
我想做的是在我的 3D 坐标系（半径 = 1）的中心定义一个球体，将圆柱形行星贴图包裹到球体表面（即在球体上执行纹理贴图），然后绘制物体周围的 3D 轨迹（如卫星轨迹）。有什么方法可以使用 matplotlib 或 mayavi 做到这一点吗？
使用 mayavi.mlab.plot3d 绘制轨迹很容易，一旦你有了自己的行星，所以我将专注于使用 mayavi 将行星映射到球体的纹理。 （原则上我们可以使用matplotlib完成任务，但是性能和质量比mayavi差很多，见本回答最后。）
精彩场景：球上球
事实证明，如果您想将球形参数化图像映射到球体上，您必须弄脏手并使用一些裸露的 vtk，但实际上需要做的工作很少，而且结果看起来不错.我将使用 Blue Marble image from NASA 进行演示。他们的自述文件说这些图片有
a geographic (Plate Carrée) projection, which is based on an equal latitude-
longitude grid spacing (not an equal area projection!)
在维基百科上查了一下，原来这也被称为equirectangular projection。也就是说，x沿线的像素直接对应经度，y沿线的像素直接对应纬度。这就是我所说的 "spherically parametrized".
所以在这种情况下，我们可以使用低级 TexturedSphereSource 来生成一个可以映射纹理的球体。自己构建球体网格可能会导致映射中出现伪影（稍后会详细介绍）。
对于底层的 vtk 工作，我重新做了 the official example found here。这就是它所需要的一切：
from mayavi import mlab
from tvtk.api import tvtk # python wrappers for the C++ vtk ecosystem
def auto_sphere(image_file):
# create a figure window (and scene)
fig = mlab.figure(size=(600, 600))
# load and map the texture
img = tvtk.JPEGReader()
img.file_name = image_file
texture = tvtk.Texture(input_connection=img.output_port, interpolate=1)
# (interpolate for a less raster appearance when zoomed in)
# use a TexturedSphereSource, a.k.a. getting our hands dirty
R = 1
Nrad = 180
# create the sphere source with a given radius and angular resolution
sphere = tvtk.TexturedSphereSource(radius=R, theta_resolution=Nrad,
phi_resolution=Nrad)
# assemble rest of the pipeline, assign texture
sphere_mapper = tvtk.PolyDataMapper(input_connection=sphere.output_port)
sphere_actor = tvtk.Actor(mapper=sphere_mapper, texture=texture)
fig.scene.add_actor(sphere_actor)
if __name__ == "__main__":
image_file = 'blue_marble_spherical.jpg'
auto_sphere(image_file)
mlab.show()
结果正是我们所期望的：
不太好的场景：不是球体
不幸的是，我不知道如何通过上述方法使用非球面贴图。此外，我们可能不想映射到一个完美的球体上，而是映射到一个椭圆体或类似的圆形物体上。对于这种情况，我们可能必须自己构建表面并尝试对其进行纹理映射。剧透警告：它不会那么漂亮。
从手动生成的球体开始，我们可以像以前一样加载纹理，并使用由 mlab.mesh:
构造的更高级别的对象
import numpy as np
from mayavi import mlab
from tvtk.api import tvtk
import matplotlib.pyplot as plt # only for manipulating the input image
def manual_sphere(image_file):
# caveat 1: flip the input image along its first axis
img = plt.imread(image_file) # shape (N,M,3), flip along first dim
outfile = image_file.replace('.jpg', '_flipped.jpg')
# flip output along first dim to get right chirality of the mapping
img = img[::-1,...]
plt.imsave(outfile, img)
image_file = outfile # work with the flipped file from now on
# parameters for the sphere
R = 1 # radius of the sphere
Nrad = 180 # points along theta and phi
phi = np.linspace(0, 2 * np.pi, Nrad) # shape (Nrad,)
theta = np.linspace(0, np.pi, Nrad) # shape (Nrad,)
phigrid,thetagrid = np.meshgrid(phi, theta) # shapes (Nrad, Nrad)
# compute actual points on the sphere
x = R * np.sin(thetagrid) * np.cos(phigrid)
y = R * np.sin(thetagrid) * np.sin(phigrid)
z = R * np.cos(thetagrid)
# create figure
mlab.figure(size=(600, 600))
# create meshed sphere
mesh = mlab.mesh(x,y,z)
mesh.actor.actor.mapper.scalar_visibility = False
mesh.actor.enable_texture = True # probably redundant assigning the texture later
# load the (flipped) image for texturing
img = tvtk.JPEGReader(file_name=image_file)
texture = tvtk.Texture(input_connection=img.output_port, interpolate=0, repeat=0)
mesh.actor.actor.texture = texture
# tell mayavi that the mapping from points to pixels happens via a sphere
mesh.actor.tcoord_generator_mode = 'sphere' # map is already given for a spherical mapping
cylinder_mapper = mesh.actor.tcoord_generator
# caveat 2: if prevent_seam is 1 (default), half the image is used to map half the sphere
cylinder_mapper.prevent_seam = 0 # use 360 degrees, might cause seam but no fake data
#cylinder_mapper.center = np.array([0,0,0]) # set non-trivial center for the mapping sphere if necessary
正如您在代码中看到的注释，有一些注意事项。第一个是球形映射模式出于某种原因翻转了输入图像（这导致反射地球）。因此，使用这种方法，我们首先必须创建输入图像的翻转版本。每张图片只需执行一次，但我在上述函数的顶部留下了相应的代码块。
第二个注意事项是，如果纹理贴图器的 prevent_seam 属性保留默认值 1，则贴图发生在 0 到 180 度方位角，球体的另一半得到反射贴图。我们显然不希望这样：我们想要绘制从 0 到 360 度方位角的整个球体。碰巧的是，此映射可能意味着我们在 phi=0 处的映射中看到接缝（不连续性），即在地图的边缘。这是尽可能使用第一种方法的另一个原因。无论如何，这是结果，包含 phi=0 点（证明没有接缝）：
圆柱映射
上述球形映射的工作方式是，表面上的每个点都通过 space 中的给定点投影到球体上。对于第一个例子，这个点是原点，对于第二个例子，我们可以设置一个长度为 3 的数组作为 cylinder_mapper.center 的值，以便映射到非原点中心的球体上。
现在，您的问题提到了圆柱映射。原则上我们可以使用第二种方法来做到这一点：
mesh.actor.tcoord_generator_mode = 'cylinder'
cylinder_mapper = mesh.actor.tcoord_generator
cylinder_mapper.automatic_cylinder_generation = 0 # use manual cylinder from points
cylinder_mapper.point1 = np.array([0,0,-R])
cylinder_mapper.point2 = np.array([0,0,R])
cylinder_mapper.prevent_seam = 0 # use 360 degrees, causes seam but no fake data
这会将球形映射更改为圆柱形。它根据设置圆柱轴和范围的两个点（[0,0,-R] 和 [0,0,R]）定义投影。每个点都根据其圆柱坐标(phi,z)映射：0到360度的方位角和坐标的垂直投影。先前关于接缝的评论仍然适用。
但是，如果我必须做这样的柱面映射，我肯定会尝试使用第一种方法。在最坏的情况下，这意味着我们必须将圆柱形参数化地图转换为球形参数化地图。同样，每张地图只需执行一次，并且可以使用 2d 插值轻松完成，例如使用 scipy.interpolate.RegularGridInterpolator。具体的变换你得知道你的非球面投影的具体情况，但是把它变换成球面投影应该不会太难，然后你就可以根据第一种情况使用TexturedSphereSource。
附录：matplotlib
为了完整起见，你可以使用 matplotlib 做你想做的事，但它会占用更多的内存和 CPU（注意你必须使用 mayavi 或 matplotlib，但你可以'不要将两者混合在一个图中）。这个想法是定义一个对应于输入地图像素的网格，并将图像作为 Axes3D.plot_surface 的 facecolors 关键字参数传递。该构造使得球体的分辨率直接耦合到映射的分辨率。我们只能使用少量的点来保持内存的易处理性，但结果看起来像素化很严重。不管怎样：
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
def mpl_sphere(image_file):
img = plt.imread(image_file)
# define a grid matching the map size, subsample along with pixels
theta = np.linspace(0, np.pi, img.shape[0])
phi = np.linspace(0, 2*np.pi, img.shape[1])
count = 180 # keep 180 points along theta and phi
theta_inds = np.linspace(0, img.shape[0] - 1, count).round().astype(int)
phi_inds = np.linspace(0, img.shape[1] - 1, count).round().astype(int)
theta = theta[theta_inds]
phi = phi[phi_inds]
img = img[np.ix_(theta_inds, phi_inds)]
theta,phi = np.meshgrid(theta, phi)
R = 1
# sphere
x = R * np.sin(theta) * np.cos(phi)
y = R * np.sin(theta) * np.sin(phi)
z = R * np.cos(theta)
# create 3d Axes
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_surface(x.T, y.T, z.T, facecolors=img/255, cstride=1, rstride=1) # we've already pruned ourselves
# make the plot more spherical
ax.axis('scaled')
if __name__ == "__main__":
image_file = 'blue_marble.jpg'
mpl_sphere(image_file)
plt.show()
上面的count参数定义了地图的下采样和渲染球体的相应大小。通过上面的180设置我们得到下图：
此外，matplotlib 使用 2d 渲染器，这意味着对于复杂的 3d 对象渲染通常会以奇怪的伪像结束（特别是，扩展对象可以完全在彼此前面或后面，因此互锁几何图形通常看起来破碎的）。考虑到这些，我肯定会使用 mayavi 来绘制带纹理的球体。 （虽然matplotlib案例中的映射是在表面上逐面工作的，所以它可以直接应用于任意表面。）
我想做的是在我的 3D 坐标系（半径 = 1）的中心定义一个球体，将圆柱形行星贴图包裹到球体表面（即在球体上执行纹理贴图），然后绘制物体周围的 3D 轨迹（如卫星轨迹）。有什么方法可以使用 matplotlib 或 mayavi 做到这一点吗？
使用 mayavi.mlab.plot3d 绘制轨迹很容易，一旦你有了自己的行星，所以我将专注于使用 mayavi 将行星映射到球体的纹理。 （原则上我们可以使用matplotlib完成任务，但是性能和质量比mayavi差很多，见本回答最后。）
精彩场景：球上球
事实证明，如果您想将球形参数化图像映射到球体上，您必须弄脏手并使用一些裸露的 vtk，但实际上需要做的工作很少，而且结果看起来不错.我将使用 Blue Marble image from NASA 进行演示。他们的自述文件说这些图片有
a geographic (Plate Carrée) projection, which is based on an equal latitude- longitude grid spacing (not an equal area projection!)
在维基百科上查了一下，原来这也被称为equirectangular projection。也就是说，x沿线的像素直接对应经度，y沿线的像素直接对应纬度。这就是我所说的 "spherically parametrized".
所以在这种情况下，我们可以使用低级 TexturedSphereSource 来生成一个可以映射纹理的球体。自己构建球体网格可能会导致映射中出现伪影（稍后会详细介绍）。
对于底层的 vtk 工作，我重新做了 the official example found here。这就是它所需要的一切：
from mayavi import mlab
from tvtk.api import tvtk # python wrappers for the C++ vtk ecosystem
def auto_sphere(image_file):
# create a figure window (and scene)
fig = mlab.figure(size=(600, 600))
# load and map the texture
img = tvtk.JPEGReader()
img.file_name = image_file
texture = tvtk.Texture(input_connection=img.output_port, interpolate=1)
# (interpolate for a less raster appearance when zoomed in)
# use a TexturedSphereSource, a.k.a. getting our hands dirty
R = 1
Nrad = 180
# create the sphere source with a given radius and angular resolution
sphere = tvtk.TexturedSphereSource(radius=R, theta_resolution=Nrad,
phi_resolution=Nrad)
# assemble rest of the pipeline, assign texture
sphere_mapper = tvtk.PolyDataMapper(input_connection=sphere.output_port)
sphere_actor = tvtk.Actor(mapper=sphere_mapper, texture=texture)
fig.scene.add_actor(sphere_actor)
if __name__ == "__main__":
image_file = 'blue_marble_spherical.jpg'
auto_sphere(image_file)
mlab.show()
结果正是我们所期望的：
不太好的场景：不是球体
不幸的是，我不知道如何通过上述方法使用非球面贴图。此外，我们可能不想映射到一个完美的球体上，而是映射到一个椭圆体或类似的圆形物体上。对于这种情况，我们可能必须自己构建表面并尝试对其进行纹理映射。剧透警告：它不会那么漂亮。
从手动生成的球体开始，我们可以像以前一样加载纹理，并使用由 mlab.mesh:
import numpy as np
from mayavi import mlab
from tvtk.api import tvtk
import matplotlib.pyplot as plt # only for manipulating the input image
def manual_sphere(image_file):
# caveat 1: flip the input image along its first axis
img = plt.imread(image_file) # shape (N,M,3), flip along first dim
outfile = image_file.replace('.jpg', '_flipped.jpg')
# flip output along first dim to get right chirality of the mapping
img = img[::-1,...]
plt.imsave(outfile, img)
image_file = outfile # work with the flipped file from now on
# parameters for the sphere
R = 1 # radius of the sphere
Nrad = 180 # points along theta and phi
phi = np.linspace(0, 2 * np.pi, Nrad) # shape (Nrad,)
theta = np.linspace(0, np.pi, Nrad) # shape (Nrad,)
phigrid,thetagrid = np.meshgrid(phi, theta) # shapes (Nrad, Nrad)
# compute actual points on the sphere
x = R * np.sin(thetagrid) * np.cos(phigrid)
y = R * np.sin(thetagrid) * np.sin(phigrid)
z = R * np.cos(thetagrid)
# create figure
mlab.figure(size=(600, 600))
# create meshed sphere
mesh = mlab.mesh(x,y,z)
mesh.actor.actor.mapper.scalar_visibility = False
mesh.actor.enable_texture = True # probably redundant assigning the texture later
# load the (flipped) image for texturing
img = tvtk.JPEGReader(file_name=image_file)
texture = tvtk.Texture(input_connection=img.output_port, interpolate=0, repeat=0)
mesh.actor.actor.texture = texture
# tell mayavi that the mapping from points to pixels happens via a sphere
mesh.actor.tcoord_generator_mode = 'sphere' # map is already given for a spherical mapping
cylinder_mapper = mesh.actor.tcoord_generator
# caveat 2: if prevent_seam is 1 (default), half the image is used to map half the sphere
cylinder_mapper.prevent_seam = 0 # use 360 degrees, might cause seam but no fake data
#cylinder_mapper.center = np.array([0,0,0]) # set non-trivial center for the mapping sphere if necessary
正如您在代码中看到的注释，有一些注意事项。第一个是球形映射模式出于某种原因翻转了输入图像（这导致反射地球）。因此，使用这种方法，我们首先必须创建输入图像的翻转版本。每张图片只需执行一次，但我在上述函数的顶部留下了相应的代码块。
第二个注意事项是，如果纹理贴图器的 prevent_seam 属性保留默认值 1，则贴图发生在 0 到 180 度方位角，球体的另一半得到反射贴图。我们显然不希望这样：我们想要绘制从 0 到 360 度方位角的整个球体。碰巧的是，此映射可能意味着我们在 phi=0 处的映射中看到接缝（不连续性），即在地图的边缘。这是尽可能使用第一种方法的另一个原因。无论如何，这是结果，包含 phi=0 点（证明没有接缝）：
圆柱映射
上述球形映射的工作方式是，表面上的每个点都通过 space 中的给定点投影到球体上。对于第一个例子，这个点是原点，对于第二个例子，我们可以设置一个长度为 3 的数组作为 cylinder_mapper.center 的值，以便映射到非原点中心的球体上。
现在，您的问题提到了圆柱映射。原则上我们可以使用第二种方法来做到这一点：
mesh.actor.tcoord_generator_mode = 'cylinder'
cylinder_mapper = mesh.actor.tcoord_generator
cylinder_mapper.automatic_cylinder_generation = 0 # use manual cylinder from points
cylinder_mapper.point1 = np.array([0,0,-R])
cylinder_mapper.point2 = np.array([0,0,R])
cylinder_mapper.prevent_seam = 0 # use 360 degrees, causes seam but no fake data
这会将球形映射更改为圆柱形。它根据设置圆柱轴和范围的两个点（[0,0,-R] 和 [0,0,R]）定义投影。每个点都根据其圆柱坐标(phi,z)映射：0到360度的方位角和坐标的垂直投影。先前关于接缝的评论仍然适用。
但是，如果我必须做这样的柱面映射，我肯定会尝试使用第一种方法。在最坏的情况下，这意味着我们必须将圆柱形参数化地图转换为球形参数化地图。同样，每张地图只需执行一次，并且可以使用 2d 插值轻松完成，例如使用 scipy.interpolate.RegularGridInterpolator。具体的变换你得知道你的非球面投影的具体情况，但是把它变换成球面投影应该不会太难，然后你就可以根据第一种情况使用TexturedSphereSource。
附录：matplotlib
为了完整起见，你可以使用 matplotlib 做你想做的事，但它会占用更多的内存和 CPU（注意你必须使用 mayavi 或 matplotlib，但你可以'不要将两者混合在一个图中）。这个想法是定义一个对应于输入地图像素的网格，并将图像作为 Axes3D.plot_surface 的 facecolors 关键字参数传递。该构造使得球体的分辨率直接耦合到映射的分辨率。我们只能使用少量的点来保持内存的易处理性，但结果看起来像素化很严重。不管怎样：
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
def mpl_sphere(image_file):
img = plt.imread(image_file)
# define a grid matching the map size, subsample along with pixels
theta = np.linspace(0, np.pi, img.shape[0])
phi = np.linspace(0, 2*np.pi, img.shape[1])
count = 180 # keep 180 points along theta and phi
theta_inds = np.linspace(0, img.shape[0] - 1, count).round().astype(int)
phi_inds = np.linspace(0, img.shape[1] - 1, count).round().astype(int)
theta = theta[theta_inds]
phi = phi[phi_inds]
img = img[np.ix_(theta_inds, phi_inds)]
theta,phi = np.meshgrid(theta, phi)
R = 1
# sphere
x = R * np.sin(theta) * np.cos(phi)
y = R * np.sin(theta) * np.sin(phi)
z = R * np.cos(theta)
# create 3d Axes
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_surface(x.T, y.T, z.T, facecolors=img/255, cstride=1, rstride=1) # we've already pruned ourselves
# make the plot more spherical
ax.axis('scaled')
if __name__ == "__main__":
image_file = 'blue_marble.jpg'
mpl_sphere(image_file)
plt.show()
上面的count参数定义了地图的下采样和渲染球体的相应大小。通过上面的180设置我们得到下图：
此外，matplotlib 使用 2d 渲染器，这意味着对于复杂的 3d 对象渲染通常会以奇怪的伪像结束（特别是，扩展对象可以完全在彼此前面或后面，因此互锁几何图形通常看起来破碎的）。考虑到这些，我肯定会使用 mayavi 来绘制带纹理的球体。 （虽然matplotlib案例中的映射是在表面上逐面工作的，所以它可以直接应用于任意表面。）