Master 3D Visualization and Mesh Processing with Python
PyVista is a powerful 3D visualization and mesh analysis library built on top of the Visualization Toolkit (VTK). It provides a high-level, Pythonic interface for scientific visualization, making complex 3D rendering and mesh processing accessible to Python developers. PyVista dramatically simplifies the VTK API while maintaining full access to its powerful capabilities.
PyVista is used across diverse scientific and engineering fields:
PyVista requires Python 3.8 or later and works on Windows, macOS, and Linux. Key dependencies include:
Installing PyVista is straightforward using pip or conda. Choose the method that works best for your environment:
# Using pip (recommended for most users)
pip install pyvista
# Using conda (for Anaconda/Miniconda users)
conda install -c conda-forge pyvista
# Install with all optional dependencies
pip install pyvista[all]
# For Jupyter notebook support only
pip install pyvista[jupyter]
# For development (includes testing tools)
pip install pyvista[dev]After installation, verify that PyVista is working correctly:
import pyvista as pv
import numpy as np
# Check version
print(f"PyVista version: {pv.__version__}")
# Test basic functionality
sphere = pv.Sphere()
print(f"Created sphere with {sphere.n_points} points")
# Check GPU support
print(f"OpenGL version: {pv.Report().opengl_version}")Configure PyVista's global settings to customize the appearance and behavior of your visualizations:
import pyvista as pv
import numpy as np
# Set default plot theme
pv.set_plot_theme('document') # Clean white background
# Other themes: 'dark', 'paraview', 'night', 'default'
# Configure global settings
pv.global_theme.font.family = 'arial'
pv.global_theme.font.size = 12
pv.global_theme.font.color = 'black'
# Window defaults
pv.global_theme.window_size = [1024, 768]
pv.global_theme.background = 'white'
# Enable anti-aliasing for smoother rendering
pv.global_theme.multi_samples = 8
# Lighting defaults
pv.global_theme.lighting = True
pv.global_theme.smooth_shading = TrueLet's create your first 3D visualization using PyVista. We'll start simple and progressively add features:
# Method 1: Quick plot (simplest)
sphere = pv.Sphere(radius=1.0)
sphere.plot()
# Method 2: With basic customization
sphere.plot(
color='lightblue',
show_edges=True,
edge_color='black',
lighting=True,
smooth_shading=True
)
# Method 3: Using Plotter for full control
plotter = pv.Plotter()
plotter.add_mesh(sphere, color='lightblue', show_edges=True)
plotter.add_text('My First PyVista Plot', font_size=14)
plotter.show_axes()
plotter.camera_position = 'iso' # Isometric view
plotter.show()PyVista provides several mesh types for different data structures:
Best for surfaces, point clouds, and polygonal data:
# Create surface mesh
sphere = pv.Sphere() # Returns PolyData
print(f"Type: {type(sphere)}") # pyvista.PolyData
# PolyData is ideal for:
# - Surface meshes (STL files)
# - Point clouds
# - Wireframe models
# - 2D and 3D polygonal dataFor 3D cells with irregular topology (FEA, CFD):
# Common in finite element analysis
# Supports tetrahedra, hexahedra, pyramids, wedges
grid = pv.UnstructuredGrid()
# Load from FEA software
# grid = pv.read('simulation_results.vtk')Regular topology with variable geometry:
# Create wavy surface
x = np.arange(-10, 10, 0.5)
y = np.arange(-10, 10, 0.5)
x, y = np.meshgrid(x, y)
z = np.sin(np.sqrt(x**2 + y**2))
grid = pv.StructuredGrid(x, y, z)
# Ideal for terrain, deformed gridsUniform 3D grids (medical imaging, voxel data):
# Create uniform grid
grid = pv.ImageData(
dimensions=(50, 50, 50),
spacing=(0.1, 0.1, 0.1),
origin=(0, 0, 0)
)
# Perfect for CT/MRI scans, voxel artCommon operations you'll use frequently:
# Create primitives
sphere = pv.Sphere()
cube = pv.Cube()
cylinder = pv.Cylinder()
# Load from file
mesh = pv.read('model.stl')
# Get information
print(f"Points: {mesh.n_points}")
print(f"Cells: {mesh.n_cells}")
print(f"Bounds: {mesh.bounds}")
# Transform
mesh.translate([1, 0, 0])
mesh.rotate_x(45)
mesh.scale(2.0)
# Visualize
mesh.plot()pv.examples.download_*() methods to access sample datasets for learning. Try pv.examples.download_bunny(), pv.examples.download_dragon(), or pv.examples.download_cow()!
PyVista can automatically detect and read numerous mesh formats. The pv.read() function is your gateway to loading external data:
# Auto-detect format from extension
mesh = pv.read('model.stl')
mesh = pv.read('data.vtk')
mesh = pv.read('scan.ply')
mesh = pv.read('model.obj')
# Read from URL
mesh = pv.read('https://example.com/model.stl')
# Download example datasets (great for learning!)
bunny = pv.examples.download_bunny()
cow = pv.examples.download_cow()
dragon = pv.examples.download_dragon()
beethoven = pv.examples.download_beethoven()
# Read multiple files
meshes = pv.read('pattern_*.vtu') # Returns MultiBlockPyVista provides built-in functions to create common geometric shapes. These are perfect for testing, prototyping, and creating simple scenes:
# Basic sphere
sphere = pv.Sphere()
# Customized sphere
sphere = pv.Sphere(
radius=2.0,
center=(1, 2, 3),
theta_resolution=50, # Longitude divisions
phi_resolution=50 # Latitude divisions
)
# Hemisphere (half sphere)
hemisphere = pv.Sphere(
start_phi=0,
end_phi=90
)# Cube (equal sides)
cube = pv.Cube(
center=(0, 0, 0),
x_length=1.0,
y_length=1.0,
z_length=1.0
)
# Box (custom dimensions)
box = pv.Box(bounds=[-1, 1, -2, 2, -0.5, 0.5])# Cylinder
cylinder = pv.Cylinder(
center=(0, 0, 0),
direction=(0, 1, 0),
radius=1.0,
height=3.0,
resolution=50,
capping=True # Cap the ends
)
# Cone
cone = pv.Cone(
center=(0, 0, 0),
direction=(0, 0, 1),
height=2.0,
radius=1.0,
resolution=50
)# Plane
plane = pv.Plane(
center=(0, 0, 0),
direction=(0, 0, 1),
i_size=10,
j_size=10,
i_resolution=10,
j_resolution=10
)
# Disk
disk = pv.Disk(
center=(0, 0, 0),
inner=0.0,
outer=1.0,
r_res=5,
c_res=20
)Create complex mathematical surfaces:
# Torus (donut shape)
torus = pv.ParametricTorus(
ringradius=1.0,
crosssectionradius=0.3
)
# Ellipsoid
ellipsoid = pv.ParametricEllipsoid(
xradius=1.0,
yradius=2.0,
zradius=0.5
)
# Mobius strip
mobius = pv.ParametricMobius()
# Klein bottle
klein = pv.ParametricKlein()Once you have a mesh, you can query its properties to understand its structure and characteristics:
# Basic information
print(f"Type: {type(mesh)}")
print(f"Number of points: {mesh.n_points}")
print(f"Number of cells: {mesh.n_cells}")
print(f"Number of arrays: {mesh.n_arrays}")
# Geometric properties
print(f"Bounds: {mesh.bounds}") # [xmin, xmax, ymin, ymax, zmin, zmax]
print(f"Center: {mesh.center}")
print(f"Length: {mesh.length}")
# For PolyData meshes
if isinstance(mesh, pv.PolyData):
print(f"Surface area: {mesh.area}")
print(f"Volume: {mesh.volume}")
# Memory usage
print(f"Memory: {mesh.memory_size / 1e6:.2f} MB")
# Data arrays
print(f"Point arrays: {mesh.point_data.keys()}")
print(f"Cell arrays: {mesh.cell_data.keys()}")Mesh data is stored as NumPy arrays, making it easy to manipulate:
# Points (vertices) - NumPy array
points = mesh.points
print(f"Points shape: {points.shape}") # (n_points, 3)
# Get coordinates
x_coords = mesh.points[:, 0]
y_coords = mesh.points[:, 1]
z_coords = mesh.points[:, 2]
# Bounds and extents
xmin, xmax, ymin, ymax, zmin, zmax = mesh.bounds
print(f"X range: [{xmin:.2f}, {xmax:.2f}]")
print(f"Y range: [{ymin:.2f}, {ymax:.2f}]")
print(f"Z range: [{zmin:.2f}, {zmax:.2f}]")PyVista offers multiple ways to visualize meshes:
# Surface rendering (default)
mesh.plot(style='surface')
# Wireframe
mesh.plot(style='wireframe', line_width=2)
# Points
mesh.plot(style='points', point_size=10)
# With edges
mesh.plot(show_edges=True, edge_color='black')
# Multiple views comparison
p = pv.Plotter(shape=(2, 2))
p.subplot(0, 0)
p.add_mesh(mesh, style='surface', color='white')
p.add_text('Surface')
p.subplot(0, 1)
p.add_mesh(mesh, style='wireframe', color='blue')
p.add_text('Wireframe')
p.subplot(1, 0)
p.add_mesh(mesh, style='points', color='red')
p.add_text('Points')
p.subplot(1, 1)
p.add_mesh(mesh, show_edges=True, color='lightblue')
p.add_text('With Edges')
p.show()import pyvista as pv
# Download sample mesh
mesh = pv.examples.download_bunny()
# Comprehensive inspection
print("="*50)
print("MESH INFORMATION")
print("="*50)
print(f"Type: {type(mesh).__name__}")
print(f"Points: {mesh.n_points:,}")
print(f"Cells: {mesh.n_cells:,}")
print(f"\nBounds:")
print(f" X: [{mesh.bounds[0]:.3f}, {mesh.bounds[1]:.3f}]")
print(f" Y: [{mesh.bounds[2]:.3f}, {mesh.bounds[3]:.3f}]")
print(f" Z: [{mesh.bounds[4]:.3f}, {mesh.bounds[5]:.3f}]")
print(f"\nGeometry:")
print(f" Area: {mesh.area:.3f}")
print(f" Volume: {mesh.volume:.3f}")
print(f" Center: {mesh.center}")
# Visualize
mesh.plot(
color='tan',
show_edges=True,
edge_color='black',
lighting=True
)cpos='xy', cpos='xz', cpos='yz' for orthographic views, or cpos='iso' for isometric view in the plot() method.
Transform meshes in 3D space using translation, rotation, and scaling operations. These operations are fundamental for positioning and orienting meshes.
# Translate (move) mesh
translated = mesh.translate([1, 2, 3], inplace=False)
# In-place modification
mesh.translate([1, 0, 0], inplace=True)
# Center mesh at origin
centered = mesh.translate(-np.array(mesh.center), inplace=False)Rotate meshes around axes or arbitrary vectors:
# Rotate around X, Y, Z axes
rotated_x = mesh.rotate_x(45, inplace=False) # 45 degrees
rotated_y = mesh.rotate_y(30, inplace=False)
rotated_z = mesh.rotate_z(60, inplace=False)
# Rotate around arbitrary vector
vector = [1, 1, 1] # Diagonal axis
angle = 45
rotated = mesh.rotate_vector(
vector,
angle,
point=mesh.center,
inplace=False
)
# Chain rotations
rotated = mesh.rotate_x(45).rotate_y(30).rotate_z(60)Resize meshes uniformly or along specific axes:
# Uniform scaling
scaled = mesh.scale(2.0, inplace=False) # 2x in all directions
# Non-uniform scaling
scaled = mesh.scale([2, 1, 0.5], inplace=False) # Different per axis
# Scale from specific point
scaled = mesh.scale(2.0, origin=mesh.center, inplace=False)# Flip along axes
flipped_x = mesh.flip_x(inplace=False)
flipped_y = mesh.flip_y(inplace=False)
flipped_z = mesh.flip_z(inplace=False)
# Reflect across plane
normal = [0, 0, 1]
point = [0, 0, 0]
reflected = mesh.reflect(normal, point=point, inplace=False)Combine meshes using CSG (Constructive Solid Geometry) operations. These are powerful for creating complex shapes from simple primitives.
# Create overlapping shapes
sphere1 = pv.Sphere(radius=1.0, center=(0, 0, 0))
sphere2 = pv.Sphere(radius=0.8, center=(0.5, 0, 0))
# CRITICAL: Clean meshes first!
sphere1 = sphere1.clean()
sphere2 = sphere2.clean()
# Union (combine)
union = sphere1.boolean_union(sphere2)
# Difference (subtract)
difference = sphere1.boolean_difference(sphere2)
# Intersection (overlap only)
intersection = sphere1.boolean_intersection(sphere2)
# Visualize all operations
p = pv.Plotter(shape=(2, 2))
p.subplot(0, 0)
p.add_mesh(sphere1, color='red', opacity=0.5)
p.add_mesh(sphere2, color='blue', opacity=0.5)
p.add_text('Original')
p.subplot(0, 1)
p.add_mesh(union, color='green')
p.add_text('Union')
p.subplot(1, 0)
p.add_mesh(difference, color='yellow')
p.add_text('Difference')
p.subplot(1, 1)
p.add_mesh(intersection, color='purple')
p.add_text('Intersection')
p.show()mesh.clean() before boolean operations! Non-manifold edges or holes will cause failures. Meshes must be closed and manifold.
Assess and visualize mesh quality using various metrics:
# Compute cell quality
quality = mesh.compute_cell_quality(
quality_measure='scaled_jacobian'
)
# Available measures:
# 'area', 'aspect_ratio', 'condition', 'diagonal',
# 'jacobian', 'max_angle', 'min_angle', 'radius_ratio',
# 'scaled_jacobian', 'shape', 'skew', 'stretch', 'volume'
# Visualize quality
quality.plot(
scalars='CellQuality',
cmap='RdYlGn',
clim=[0, 1],
show_scalar_bar=True
)
# Find poor quality cells
bad_cells = quality.threshold(
value=0.3,
scalars='CellQuality',
invert=True
)
print(f"Poor quality cells: {bad_cells.n_cells}")Compute surface normals for lighting and analysis:
# Compute normals
mesh_with_normals = mesh.compute_normals(
cell_normals=True,
point_normals=True,
split_vertices=False,
flip_normals=False,
consistent_normals=True,
feature_angle=30
)
# Access normals
point_normals = mesh_with_normals['Normals']
# Visualize normals with arrows
arrows = mesh.glyph(
orient='Normals',
scale=False,
factor=0.1
)
p = pv.Plotter()
p.add_mesh(mesh, color='white')
p.add_mesh(arrows, color='blue')
p.show()Analyze surface curvature to identify features:
# Mean curvature
curv_mean = mesh.compute_curvature(curvature_type='mean')
# Gaussian curvature
curv_gaussian = mesh.compute_curvature(curvature_type='gaussian')
# Maximum/minimum curvature
curv_max = mesh.compute_curvature(curvature_type='maximum')
curv_min = mesh.compute_curvature(curvature_type='minimum')
# Compare all types
p = pv.Plotter(shape=(2, 2))
p.subplot(0, 0)
p.add_mesh(curv_mean, scalars='mean_curvature', cmap='coolwarm')
p.add_text('Mean Curvature')
p.subplot(0, 1)
p.add_mesh(curv_gaussian, scalars='gaussian_curvature', cmap='coolwarm')
p.add_text('Gaussian Curvature')
p.subplot(1, 0)
p.add_mesh(curv_max, scalars='maximum_curvature', cmap='coolwarm')
p.add_text('Maximum')
p.subplot(1, 1)
p.add_mesh(curv_min, scalars='minimum_curvature', cmap='coolwarm')
p.add_text('Minimum')
p.show()Fix common mesh problems and prepare meshes for analysis:
# Remove duplicates and merge points
cleaned = mesh.clean(
point_merging=True,
tolerance=1e-6,
remove_unused_points=True
)
# Fill holes
filled = mesh.fill_holes(hole_size=10.0)
# Extract surface (remove internal cells)
surface = mesh.extract_surface()
# Triangulate (convert all cells to triangles)
triangulated = mesh.triangulate()def repair_mesh(mesh):
"""Complete mesh repair pipeline."""
print("Repairing mesh...")
print(f"Initial: {mesh.n_points} points, {mesh.n_cells} cells")
# Step 1: Clean
mesh = mesh.clean()
# Step 2: Extract surface
mesh = mesh.extract_surface()
# Step 3: Fill holes
mesh = mesh.fill_holes(100.0)
# Step 4: Triangulate
mesh = mesh.triangulate()
# Step 5: Keep largest region
mesh = mesh.connectivity(largest=True)
# Step 6: Compute normals
mesh = mesh.compute_normals()
print(f"Final: {mesh.n_points} points, {mesh.n_cells} cells")
return mesh
# Use it
repaired = repair_mesh(damaged_mesh)Analyze and extract connected regions:
# Extract all connected regions
connected = mesh.connectivity('all')
# Extract largest region only
largest = mesh.connectivity(largest=True)
# Count regions
n_regions = len(np.unique(connected['RegionId']))
print(f"Connected regions: {n_regions}")
# Split into separate meshes
split_meshes = connected.split_bodies()
print(f"Split into {len(split_meshes)} separate meshes")# Add scalar data
mesh['temperature'] = np.random.rand(mesh.n_points) * 100
# Visualize
mesh.plot(
scalars='temperature',
cmap='coolwarm',
show_scalar_bar=True
)# Add vectors
vectors = np.random.rand(mesh.n_points, 3)
mesh['vectors'] = vectors
# Visualize with arrows
arrows = mesh.glyph(
orient='vectors',
scale='vectors',
factor=0.1
)
arrows.plot()# Create contours
contours = mesh.contour(
isosurfaces=5,
scalars='temperature'
)
# Slice mesh
slice_x = mesh.slice(normal='x')
clipped = mesh.clip(normal='z', value=0)# Point data
mesh['point_temp'] = np.random.rand(mesh.n_points) * 100
# Cell data
mesh['cell_quality'] = np.random.rand(mesh.n_cells)
# Vector data
velocities = np.random.rand(mesh.n_points, 3)
mesh['velocity'] = velocities# Cell to point
mesh_p = mesh.cell_data_to_point_data()
# Point to cell
mesh_c = mesh.point_data_to_cell_data()
# Interpolation
target = pv.Sphere()
interpolated = target.interpolate(mesh)# Filter by value
high_values = mesh.threshold(
value=70.0,
scalars='temperature'
)
# Range selection
mid_range = mesh.threshold(
value=[40, 60],
scalars='temperature'
)# Points as NumPy array
points = mesh.points
print(f"Shape: {points.shape}")
# Modify coordinates
mesh.points[:, 2] *= 2 # Scale Z
# Statistical operations
temp = mesh['temperature']
mean_temp = np.mean(temp)
std_temp = np.std(temp)
# Vector magnitude
v = mesh['velocity']
magnitude = np.linalg.norm(v, axis=1)
mesh['speed'] = magnitudeimport pandas as pd
# Create DataFrame
df = pd.DataFrame({
'x': mesh.points[:, 0],
'y': mesh.points[:, 1],
'z': mesh.points[:, 2],
'temp': mesh['temperature']
})
# Analysis
summary = df.describe()
correlation = df.corr()
# Filter
high_temp = df[df['temp'] > 70]
# Back to mesh
points = df[['x', 'y', 'z']].values
cloud = pv.PolyData(points)
cloud['temp'] = df['temp'].values# Load volume
volume = pv.examples.download_head()
# Render
p = pv.Plotter()
p.add_volume(
volume,
cmap='bone',
opacity='sigmoid'
)
p.show()# Create vector field
mesh['vectors'] = np.random.rand(mesh.n_points, 3)
# Generate streamlines
streamlines = mesh.streamlines(
vectors='vectors',
start_position=(0, 0, 0),
max_time=1.0
)
# Visualize
p = pv.Plotter()
p.add_mesh(mesh, opacity=0.2)
p.add_mesh(streamlines, line_width=2)
p.show()# Laplacian smoothing
smoothed = mesh.smooth(n_iter=100)
# Windowed sinc
smoothed = mesh.smooth_taubin(n_iter=20)# Reduce complexity
decimated = mesh.decimate(target_reduction=0.75)
# Preserve topology
decimated = mesh.decimate_pro(
target_reduction=0.5,
preserve_topology=True
)# Linear subdivision
subdivided = mesh.subdivide(nsub=2, subfilter='linear')
# Loop subdivision
subdivided = mesh.subdivide(nsub=2, subfilter='loop')# Create terrain
x = np.arange(0, 100, 1)
y = np.arange(0, 100, 1)
x, y = np.meshgrid(x, y)
z = 10 * np.sin(x/10) * np.cos(y/10) + 50
terrain = pv.StructuredGrid(x, y, z)
terrain['elevation'] = z.ravel()
terrain.plot(
scalars='elevation',
cmap='gist_earth',
lighting=True
)# Load scan
volume = pv.read('brain_scan.nii')
# Extract tissue
skull = volume.threshold(value=500)
brain = volume.threshold([200, 500])
# Visualize
p = pv.Plotter()
p.add_mesh(skull, color='white', opacity=0.8)
p.add_mesh(brain, color='pink', opacity=0.6)
p.show()# Slider widget
p = pv.Plotter()
p.add_mesh(mesh)
def update_opacity(value):
mesh.opacity = value
p.add_slider_widget(
update_opacity,
[0, 1],
value=0.5,
title='Opacity'
)
p.show()# Enable picking
p = pv.Plotter()
p.add_mesh(mesh, pickable=True)
def callback(point):
print(f"Picked: {point}")
p.enable_point_picking(callback)
p.show()| Format | Extension | Description |
|---|---|---|
| VTK | .vtk, .vtu, .vtp | VTK native formats |
| STL | .stl | Stereolithography |
| PLY | .ply | Polygon File Format |
| OBJ | .obj | Wavefront OBJ |
# Read
mesh = pv.read('model.stl')
# Write
mesh.save('output.vtk')
mesh.save('output.stl')
# Screenshot
p = pv.Plotter(off_screen=True)
p.add_mesh(mesh)
p.show(screenshot='output.png')# Creates copy
rotated = mesh.rotate_z(45, inplace=False)
# Modifies original
mesh.rotate_z(45) # inplace=True by default# Clean before boolean ops
mesh1 = mesh1.clean()
mesh2 = mesh2.clean()
result = mesh1.boolean_union(mesh2)decimate() for large meshesoff_screen=True for batch processing# Simplify before plotting
simplified = large_mesh.decimate(0.9)
simplified.plot()
# Off-screen rendering
p = pv.Plotter(off_screen=True)
p.add_mesh(mesh)
p.show(screenshot='output.png')p = pv.Plotter()
p.open_gif('animation.gif', fps=30)
sphere = pv.Sphere()
for angle in np.linspace(0, 360, 100):
p.clear()
rotated = sphere.rotate_z(angle, inplace=False)
p.add_mesh(rotated, color='lightblue')
p.write_frame()
p.close()p = pv.Plotter()
p.add_mesh(mesh)
path = p.generate_orbital_path(n_points=100)
p.open_gif('orbit.gif')
p.orbit_on_path(path, write_frames=True)
p.close()# Set backend
pv.set_jupyter_backend('pythreejs') # Interactive
# or
pv.set_jupyter_backend('static') # Static images
# Plot in notebook
mesh.plot()static for quick previewspythreejs for interactionipyvtklink for large datasets# Reusable function
def visualize_field(mesh, field_name, cmap='viridis'):
"""Visualize scalar field on mesh."""
p = pv.Plotter()
p.add_mesh(
mesh,
scalars=field_name,
cmap=cmap,
show_scalar_bar=True
)
return p
# Use it
plotter = visualize_field(my_mesh, 'temperature', 'plasma')
plotter.show()# Always handle file I/O
try:
mesh = pv.read('model.stl')
except FileNotFoundError:
print("File not found!")
mesh = pv.Sphere()
# Validate mesh
if mesh.n_points == 0:
raise ValueError("Empty mesh!")
# Check for NaN
if np.isnan(mesh.points).any():
mesh.points = np.nan_to_num(mesh.points)You now have comprehensive knowledge of PyVista for 3D visualization and mesh processing!