python如何表示高维数据 python高维数据可视化

python 科学计算三维可视化笔记 第三周 高级进阶

• python 科学计算三维可视化笔记 第三周 高级进阶

• 一、Mayavi 入门

• 1. Mayavi 库的基本元素
• 2. 快速绘制实例
• 3. Mayavi 管线

• 二、Mlab 基础

• （一）基于 numpy 数组的绘图函数

• 1. 0D 数据：``points3d()`` 函数
• 2. 1D 数据：``plot3d()`` 函数
• 3. 2D 数据：``imshow()`` 函数
• 4. 2D 数据：``surf()`` 函数
• 5. 2D 数据：``contour_surf()`` 函数
• 6. 2D 数据：其他函数
• 7. 3D 数据：``contour3d()`` 函数
• 8. 3D 数据：``quiver3d()`` 函数
• 9. 3D 数据：``flow()`` 函数

• （二）改变物体的颜色
• （三）mlab 控制函数
• （四）鼠标选取交互操作

• 1. ``on_mouse_pick()`` 函数
• 2. 选取红色小球问题分析
• 3. 代码
• 4. 程序框架的优化
• 5. 优化代码
• 6. 结果

• （五）mlab 管线控制函数
• （六）实例1：标量数据可视化

• 1. 绘制等值面
• 2. 绘制切平面
• 3. 复合观测方法

• （七）实例2：矢量数据可视化

• 1. 使用 ``quiver3d()`` 函数
• 2. masking vector 采样
• 3. Cut Plane 切面
• 4. 级数的等值面
• 5. 流线 flow 的可视化
• 6. 复合观测方法

---

内容来自中国大学MOOC，北京理工大学，python数据分析与展示课程，侵删。
如有错误，烦请指出。

---

python 科学计算三维可视化笔记 第三周 高级进阶

一、Mayavi 入门

1. Mayavi 库的基本元素

Mayavi.mlab 模块：

Mayavi API：

2. 快速绘制实例

(1) 实例1

'''快速绘制实例1'''
from mayavi import mlab

# 建立数据
x = [[-1,1,1,-1,-1], [-1,1,1,-1,-1]]
y = [[-1,-1,-1,-1,-1], [1,1,1,1,1]]
z = [[1,1,-1,-1,1], [1,1,-1,-1,1]]

# 对数据进行三维可视化
s = mlab.mesh(x,y,z)
mlab.show()

  

(2) 实例2

'''快速绘制实例2'''
from numpy import pi, sin, cos, mgrid
from mayavi import mlab

# 建立数据
dphi, dtheta = pi/250.0, pi/250.0
[phi,theta] = mgrid[0:pi+dphi*1.5:dphi,0:2*pi+dtheta*1.5:dtheta]
m0 = 4; m1 = 3; m2 = 2; m3 = 3; m4 = 6; m5 = 2; m6 = 6; m7 = 4;
r = sin(m0*phi)**m1 + cos(m2*phi)**m3 + sin(m4*theta)**m5 + cos(m6*theta)**m7
x = r*sin(phi)*cos(theta)
y = r*cos(phi)
z = r*sin(phi)*sin(theta)

# 对数据进行三维可视化
s = mlab.mesh(x, y, z)
# s = mlab.mesh(x, y, z, representation='wireframe', line_width=1.0)
mlab.show()

mlab.mesh() 函数默认为 surface（左图），可以更改为 wireframe（右图）：

  

可以用键盘和鼠标对场景进行操作：

• 旋转场景：左键拖动或键盘的方向键
• 平移场景：按住 shift 键并使用左键拖动，shift+方向键盘
• 缩放场景：鼠标右键上下拖动或使用 + 和 - 按键
• 滚动相机：按住 ctrl 键并用左键拖动
• 工具栏：从坐标轴 6 个方向观察场景、等角投影、切换平行透视和成角透视等

3. Mayavi 管线

• engine：建立和销毁 scenes
• scenes：多个数据集合 sources
• filters：对数据进行变换
• module manager：控制颜色，Colors and legends
• modules：最终数据的表示，如线条、平面等

点击 Mayavi pipeline 或使用 mlab.show_pipeline() 打开 mayavi 管线 对话框：

(1) 管线中的对象 scene

Mayavi Scene：处于树的最顶层的对象，表示场景

使用 s = mlab.gcf() 获取当前场景

使用 print(s.scene.background) 输出当前场景的背景色

(2) 管线中的对象 GridSource

GridSource 是网格数据源，每一项是标量数据的名称

使用 source = s.children[0] 获取 GridSource 对象

使用 print(source.name) 输出名称

使用 print(repr(source.data.points)) 输出坐标

使用 print(repr(source.data.point_data.scalars)) 输出每个点对应的标量数组

(3) 管线中的对象 PolyDataNormals

PolyDataNormals：数据源的法向量

(4) 管线中的对象 Colors and legends更改 lut mode 颜色为蓝色，勾选 show legend 添加颜色条显示颜色和标量之间的关系：

colors = manager.children[0]
colors.scalar_lut_manager.lut_mode = 'Blues'
colors.scalar_lut_manager.show_legend = 'True'

结果：

(5) 管线中的对象 Surface

将 GridSource 中输出的 PolyData 数据转换为最终在场景中显示的三维实体，将 representation 改为 wireframe，将 line width 设置为 0，将 opacity (透明度) 设置为 0.6：

surface = colors.children[0]
surface.actor.property.representation = 'wireframe'
surface.actor.property.opacity = 0.6

结果：

程序配置属性的步骤：

• 获取场景对象，mlab.gcf()
• 通过 chirldren 属性，在管线中找到需要修改的对象
• 配置窗口有多个选项卡，属性需要一级一级获得
• 通过工具栏 “红色原点” 菜单，可以记录配置过程

二、Mlab 基础

（一）基于 numpy 数组的绘图函数

mlab 对 numpy 建立可视化过程：

• 建立数据源
• 使用 Filter（可选）
• 添加可视化模块

3D 绘图函数，包括 points3d() 和 plot3d()：

1. 0D 数据：points3d() 函数

(1) 函数形式

points3d(x, y, z ...)
points3d(x, y, z, s, ...)
points3d(x, y, z, f, ...)

• x, y, z 表示 numpy 数组、列表或者其他形式的点三维坐标
• s 表示在该坐标点处的标量值
• f 表示通过函数 f(x, y, z) 返回的标量值

(2) 函数参数

(3) 实例

'''points3d()'''
import numpy as np
from mayavi import mlab

# 建立数据
t = np.linspace(0, 4*np.pi, 20)
x = np.sin(2*t)
y = np.cos(t)
z = np.cos(2*t)
s = 2 + np.sin(t)

# 对数据进行可视化
points = mlab.points3d(x, y, z, s, colormap='Reds', scale_factor=.25)
# points = mlab.points3d(x, y, z, s, colormap='Greens', scale_factor=.25)
mlab.show()

  

2. 1D 数据：plot3d() 函数

(1) 函数形式

plod3d(x, y, z ...)
plot3d(x, y, z, s, ...)

• x, y, z 表示 numpy 数组、列表或者其他形式的点三维坐标
• s 表示在该坐标点处的标量值

(2) 函数参数

(3) 实例

'''plot3d'''
import numpy as np
from mayavi import mlab

# 建立数据
n_mer, n_long = 6, 11
dphi = np.pi / 1000.0
phi = np.arange(0.0, 2*np.pi+0.5*dphi, dphi)
mu = phi * n_mer
x = np.cos(mu) * (1+np.cos(n_long*mu/n_mer)*0.5)
y = np.sin(mu) * (1+np.cos(n_long*mu/n_mer)*0.5)
z = np.sin(n_long*mu/n_mer) * 0.5

# 对数据进行可视化
l = mlab.plot3d(x, y, z, np.sin(mu), tube_radius=0.025, colormap='Spectral')
mlab.show()

3. 2D 数据：imshow() 函数

将二维数组可视化为一张图像
(1) 函数形式

imshow(s, ...)

• s 是一个二维数组，s 的值使用 colormap 被映射为颜色

(2) 函数参数

(3) 实例

'''imshow()'''
import numpy as np
from mayavi import mlab

# 建立数据
s = np.random.random((10, 10))

# 对数据进行可视化
img = mlab.imshow(s, colormap='gist_earth')
mlab.show()

4. 2D 数据：surf() 函数

将二维数组可视化为一个平面，z 轴描述了数组点的高度
(1) 函数形式

surf(s, ...)
surf(x, y, s, ...)
surf(x, y, f, ...)

• x, y 是一维或二维数组，由 np.ogrid 或 np.mgrid 得到
• s 是一个高程矩阵，用二维数组表示

(2) 实例

'''surf()'''
import numpy as np
from mayavi import mlab

# 建立数据
def f(x,y):
    return np.sin(x-y) + np.cos(x+y)
x,y = np.mgrid[-7.:7.05:0.1, -5:5.05:0.05]

# 对数据进行可视化
s = mlab.surf(x, y, f)
mlab.show()

5. 2D 数据：contour_surf() 函数

将二维数组可视化为等高线，高度值由数组点的值来确定
(1) 函数形式

contour_surf(s, ...)
contour_surf(x, y, s, ...)
contour_surf(x, y, f, ...)

• x, y 是一维或二维数组，由 np.ogrid 或 np.mgrid 得到
• s 是一个高程矩阵，用二维数组表示

(2) 实例

'''contour_surf()'''
import numpy as np
from mayavi import mlab

# 建立数据
def f(x,y):
    return np.sin(x-y) + np.cos(x+y)
x,y = np.mgrid[-7.:7.05:0.1, -5:5.05:0.05]

# 对数据进行可视化
s = mlab.contour_surf(x, y, f)
mlab.show()

6. 2D 数据：其他函数

• mesh() 函数绘制由三个二维数组 x、y、z 描述坐标点的网格平面
• barchart() 函数根据二维、三维或者点云数据绘制的三位柱状图
• triangular_mesh() 函数绘制由 x、y、z 坐标点描述的三角网格面

7. 3D 数据：contour3d() 函数

三维数组定义的体数据的等值面可视化
(1) 函数形式

contour3d(scalars, ...)
contour3d(x, y, z, scalars, ...)

• scalars 网格上的数据，用三维 numpy 数组表示
• x、y、z 是三维空间坐标

(2) 函数参数

• contours：定义等值面的数量
• transparent：True 表示透明

(3) 实例

'''contour3d()'''
import numpy as np
from mayavi import mlab

# 建立数据
x, y, z = np.ogrid[-5:5:64j, -5:5:64j, -5:5:64j]
scalars = x*x + y*y + z*z

# 对数据进行可视化
obj = mlab.contour3d(scalars, contours=8, transparent=True)
mlab.show()

8. 3D 数据：quiver3d() 函数

三维矢量数组的可视化，箭头表示在该点的矢量数据
(1) 函数形式

quiver3d(u, v, w, ...)
quiver3d(x, y, z, u, v, w, ...)
quiver3d(x, y, z, f, ...)

• u、v、w 是用 numpy 数组表示的向量
• x、y、z 表示箭头的位置
• f 需要返回在给定位置 (x, y, z) 的 (u, v, w) 矢量

(2) 实例

'''quiver3d()'''
import numpy as np
from mayavi import mlab

# 建立数据
x,y,z = np.mgrid[-2:3, -2:3, -2:3]
r = np.sqrt(x**2 + y**2 + z**4)
u = y*np.sin(r) / (r+0.001)
v = -x*np.sin(r) / (r+0.001)
w = np.zeros_like(z)

# 对数据进行可视化
obj = mlab.quiver3d(x, y, z, u, v, w, line_width=3, scale_factor=1)
mlab.show()

9. 3D 数据：flow() 函数

绘制三维数组描述的向量场的粒子轨迹

（二）改变物体的颜色

• colormap 定义的颜色，也叫 LUT (Look Up Table)
• LUT 是一个 255×4 的数组，列向量表示RGBA，每个值的范围从 0 到 255

'''colormaps'''
import numpy as np
from mayavi import mlab

# 建立数据
x, y = np.mgrid[-10:10:200j, -10:10:200j]
z = 100 * np.sin(x*y) / (x*y)

# 对数据进行可视化
mlab.figure(bgcolor=(1,1,1))
surf = mlab.surf(z, colormap='cool')
mlab.show()

改变 colormap（上面代码，左图）；改变 LUT（下面代码，右图）：

  

'''LUT'''
import numpy as np
from mayavi import mlab

# 建立数据
x, y = np.mgrid[-10:10:200j, -10:10:200j]
z = 100 * np.sin(x*y) / (x*y)

# 对数据进行可视化
mlab.figure(bgcolor=(1,1,1))
surf = mlab.surf(z, colormap='cool')

# 访问surf对象的LUT
lut = surf.module_manager.scalar_lut_manager.lut.table.to_array()
# 增加透明梯度，修改alpha通道
lut[:, -1] = np.linspace(0, 255, 256)
surf.module_manager.scalar_lut_manager.lut.table = lut

mlab.show()

（三）mlab 控制函数

• 图像控制函数
• 图像装饰函数
• 相机控制函数
• 其他控制函数

（四）鼠标选取交互操作

• 选取一个物体，查看数据
• 选取物体上一点，查看数据

1. on_mouse_pick() 函数

on_mouse_pick(callback, type='point', Button='Left', Remove=False)

• type：'point'、'cell'、'world'
• Button：'left'、'Middle'、'Right'
• Remove：如果值为 True，则 callback 函数不起作用
• 返回：一个 vtk picker 对象

2. 选取红色小球问题分析

• 建立一个 figure
• 随机生成红、白小球
• 初始化红色小球选取外框
• 鼠标选取任意红色小球，外框移动到该小球上 (callback)
• 建立 on_mouse_pick() 响应机制

3. 代码

'''鼠标选取'''
import numpy as np
from mayavi import mlab

####### 场景初始化 #######
figure = mlab.gcf()

# 用mlab.points3d建立红色小球的集合
x1, y1, z1 = np.random.random((3, 10))
red_glyphs = mlab.points3d(x1, y1, z1, color=(1, 0, 0), resolution=10)

# 用mlab.points3d建立白色小球的集合
x2, y2, z2 = np.random.random((3, 10))
white_glyphs = mlab.points3d(x2, y2, z2, color=(0.9, 0.9, 0.9), resolution=10)

# 绘制选取框，并放在第一个小球上
outline = mlab.outline(line_width=3)
outline.outline_mode = 'cornered'
outline.bounds = (x1[0] - 0.1, x1[0] + 0.1,
                  y1[0] - 0.1, y1[0] + 0.1,
                  z1[0] - 0.1, z1[0] + 0.1)

####### 处理选取事件 #######

# 获取构成一个红色小球的顶点列表
glyph_points = red_glyphs.glyph.glyph_source.glyph_source.output.points.to_array()

# 处理选取时间
def picker_callback(picker):
    if picker.actor in red_glyphs.actor.actors:
        # 计算被选取的小球的ID号
        point_id = int(picker.point_id / glyph_points.shape[0])  # int向下取整        
        if point_id != -1:    # 如果没有小球被选取，则point_id = -1
            # 计算与此红色小球相关的坐标
            x, y, z = x1[point_id], y1[point_id], z1[point_id]
            # 将外框移到小球上
            outline.bounds = (x - 0.1, x + 0.1,
                              y - 0.1, y + 0.1,
                              z - 0.1, z + 0.1)

# 建立响应机制
picker = figure.on_mouse_pick(picker_callback)

# 设置窗口的标题并显示窗口
mlab.title('Click on red balls')
mlab.show()

4. 程序框架的优化

• 程序运行有两个问题：小球初始速度太慢、鼠标选取不精确
• 加入以下三条语句：

figure.scene.disable_render = True
figure.scene.disable_render = False
Picker.tolerance = 0.01

5. 优化代码

'''鼠标选取 优化'''
import numpy as np
from mayavi import mlab

####### 场景初始化 #######
figure = mlab.gcf()

figure.scene.disable_render = True

# 用mlab.points3d建立红色小球的集合
x1, y1, z1 = np.random.random((3, 10))
red_glyphs = mlab.points3d(x1, y1, z1, color=(1, 0, 0), resolution=10)

# 用mlab.points3d建立白色小球的集合
x2, y2, z2 = np.random.random((3, 10))
white_glyphs = mlab.points3d(x2, y2, z2, color=(0.9, 0.9, 0.9), resolution=10)

figure.scene.disable_render = False

# 绘制选取框，并放在第一个小球上
outline = mlab.outline(line_width=3)
outline.outline_mode = 'cornered'
outline.bounds = (x1[0] - 0.1, x1[0] + 0.1,
                  y1[0] - 0.1, y1[0] + 0.1,
                  z1[0] - 0.1, z1[0] + 0.1)

####### 处理选取事件 #######

# 获取构成一个红色小球的顶点列表
glyph_points = red_glyphs.glyph.glyph_source.glyph_source.output.points.to_array()

# 处理选取时间
def picker_callback(picker):
    if picker.actor in red_glyphs.actor.actors:
        # 计算被选取的小球的ID号
        point_id = int(picker.point_id / glyph_points.shape[0])  # int向下取整        
        if point_id != -1:    # 如果没有小球被选取，则point_id = -1
            # 计算与此红色小球相关的坐标
            x, y, z = x1[point_id], y1[point_id], z1[point_id]
            # 将外框移到小球上
            outline.bounds = (x - 0.1, x + 0.1,
                              y - 0.1, y + 0.1,
                              z - 0.1, z + 0.1)

# 建立响应机制
picker = figure.on_mouse_pick(picker_callback)

# 设置 tolerance 参数，提高选取精度
Picker.tolerance = 0.01

# 设置窗口的标题并显示窗口
mlab.title('Click on red balls')
mlab.show()

6. 结果

（五）mlab 管线控制函数

mlab.pipeline.function()

• Sources：数据源
• Filters：用来数据变换
• Modules：用来实现可视化

Sources：

Filters：

Modules：

（六）实例1：标量数据可视化

1. 绘制等值面

使用 mlab.contour3d(s) 绘制等值面，只能看到最外层的等值面（左图）

'''实例1：标量数据可视化 等值面绘制'''
import numpy as np
from mayavi import mlab

# 生成标量数据
x, y, z = np.ogrid[-10:10:20j, -10:10:20j, -10:10:20j]
s = np.sin(x*y*z)/(x*y*z)

# 对数据进行可视化
mlab.contour3d(s)
mlab.show()

2. 绘制切平面

'''实例1：标量数据可视化 切平面绘制'''
from mayavi import mlab
from mayavi.tools import pipeline

mlab.pipeline.image_plane_widget(mlab.pipeline.scalar_field(s),
                                plane_orientation='x_axes',
                                slice_index=10,)
mlab.pipeline.image_plane_widget(mlab.pipeline.scalar_field(s),
                                plane_orientation='y_axes',
                                slice_index=10,)
mlab.outline()
mlab.show()

3. 复合观测方法

'''实例1：标量数据可视化 复合观测方法'''
from mayavi import mlab
from mayavi.tools import pipeline

src = mlab.pipeline.scalar_field(s)
mlab.pipeline.iso_surface(src, contours=[s.min()+0.1*s.ptp(), ], opacity=0.1)
mlab.pipeline.iso_surface(src, contours=[s.max()-0.1*s.ptp(), ])
mlab.pipeline.image_plane_widget(src,
                            plane_orientation='z_axes',
                            slice_index=10,
                        )
mlab.show()

（七）实例2：矢量数据可视化

1. 使用 quiver3d() 函数

简单，但结果非常密集

'''实例2：矢量数据可视化 quiver3d()'''
import numpy as np
from mayavi import mlab

x, y, z = np.mgrid[0:1:20j, 0:1:20j, 0:1:20j]

u =    np.sin(np.pi*x) * np.cos(np.pi*z)
v = -2*np.sin(np.pi*y) * np.cos(2*np.pi*z)
w = np.cos(np.pi*x)*np.sin(np.pi*z) + np.cos(np.pi*y)*np.sin(2*np.pi*z)


mlab.quiver3d(u,v,w)
mlab.outline()

mlab.show()

2. masking vector 采样

降低采样密度

'''实例2：矢量数据可视化 masking vector'''
import numpy as np
from mayavi import mlab

x, y, z = np.mgrid[0:1:20j, 0:1:20j, 0:1:20j]

u =    np.sin(np.pi*x) * np.cos(np.pi*z)
v = -2*np.sin(np.pi*y) * np.cos(2*np.pi*z)
w = np.cos(np.pi*x)*np.sin(np.pi*z) + np.cos(np.pi*y)*np.sin(2*np.pi*z)

src = mlab.pipeline.vector_field(u, v, w)

mlab.pipeline.vectors(src, mask_points=10, scale_factor=2.0)

mlab.show()

3. Cut Plane 切面

'''实例2：矢量数据可视化 cut plane'''
import numpy as np
from mayavi import mlab

x, y, z = np.mgrid[0:1:20j, 0:1:20j, 0:1:20j]

u =    np.sin(np.pi*x) * np.cos(np.pi*z)
v = -2*np.sin(np.pi*y) * np.cos(2*np.pi*z)
w = np.cos(np.pi*x)*np.sin(np.pi*z) + np.cos(np.pi*y)*np.sin(2*np.pi*z)

src = mlab.pipeline.vector_field(u, v, w)

mlab.pipeline.vector_cut_plane(src, mask_points=10, scale_factor=2.0)

mlab.show()

4. 级数的等值面

显示矢量的法向等值面

'''实例2：矢量数据可视化 级数的等值面'''
import numpy as np
from mayavi import mlab

x, y, z = np.mgrid[0:1:20j, 0:1:20j, 0:1:20j]

u =    np.sin(np.pi*x) * np.cos(np.pi*z)
v = -2*np.sin(np.pi*y) * np.cos(2*np.pi*z)
w = np.cos(np.pi*x)*np.sin(np.pi*z) + np.cos(np.pi*y)*np.sin(2*np.pi*z)

src = mlab.pipeline.vector_field(u, v, w)
magnitude = mlab.pipeline.extract_vector_norm(src)

mlab.pipeline.iso_surface(magnitude, contours=[2.0,0.5])
mlab.outline()

mlab.show()

5. 流线 flow 的可视化

'''实例2：矢量数据可视化 流线flow的可视化'''
import numpy as np
from mayavi import mlab

x, y, z = np.mgrid[0:1:20j, 0:1:20j, 0:1:20j]

u =    np.sin(np.pi*x) * np.cos(np.pi*z)
v = -2*np.sin(np.pi*y) * np.cos(2*np.pi*z)
w = np.cos(np.pi*x)*np.sin(np.pi*z) + np.cos(np.pi*y)*np.sin(2*np.pi*z)

flow = mlab.flow(u, v, w, seed_scale=1, seed_resolution=5, integration_direction='both')
mlab.outline()

mlab.show()

6. 复合观测方法

'''实例2：矢量数据可视化 复合观测方法'''
import numpy as np
from mayavi import mlab

x, y, z = np.mgrid[0:1:20j, 0:1:20j, 0:1:20j]

u =    np.sin(np.pi*x) * np.cos(np.pi*z)
v = -2*np.sin(np.pi*y) * np.cos(2*np.pi*z)
w = np.cos(np.pi*x)*np.sin(np.pi*z) + np.cos(np.pi*y)*np.sin(2*np.pi*z)

src = mlab.pipeline.vector_field(u, v, w)
magnitude = mlab.pipeline.extract_vector_norm(src)

iso = mlab.pipeline.iso_surface(magnitude, contours=[2.0, ], opacity=0.3)
vec = mlab.pipeline.vectors(magnitude, mask_points=40, line_width=1,
                            color=(0.8,0.8,0.8),
                            scale_factor=4.)
flow = mlab.pipeline.streamline(magnitude, seedtype='plane',
                                seed_visible=False,
                                seed_scale=0.5,
                                seed_resolution=1,
                                linetype='ribbon')
vcp = mlab.pipeline.vector_cut_plane(magnitude, mask_points=2, 
                                     scale_factor=4,
                                     colormap='jet',
                                     plane_orientation='x_axes')
mlab.outline()
mlab.show()