CODEWITHSUNDEEP
python3dpyvista
Stéphane Laurent
Stéphane Laurent

Reputation: 84659

How to clip a PyVista mesh to a ball?

I have this PyVista isosurface:

from math import pi, cos, sin
import pyvista as pv
import numpy as np

def strange_surface(x, y, z, A, B):
    return (
        z**4 * B**2
        + 4 * x * y**2 * A * B**2
        + x * z**2 * A * B**2
        - 2 * z**4 * A
        - 4 * x * y**2 * B**2
        - x * z**2 * B**2
        + 3 * z**2 * A * B**2
        - 2 * z**4
        - x * A * B**2
        - 2 * z**2 * A
        + x * B**2
        + A * B**2
        + 2 * z**2
        - B**2
    )

# generate data grid for computing the values
X, Y, Z = np.mgrid[(-3.05):3.05:250j, (-3.05):3.05:250j, (-3.05):3.05:250j]
# create a structured grid
grid = pv.StructuredGrid(X, Y, Z)
# compute and assign the values
A = cos(2.5*pi/4)
B = sin(2.5*pi/4)
values = strange_surface(X, Y, Z, A, B)
grid.point_data["values"] = values.ravel(order="F")
isosurf = grid.contour(isosurfaces=[0])
#
mesh = isosurf.extract_geometry()
dists = np.linalg.norm(mesh.points, axis=1)
dists = (dists - dists.min()) / (dists.max() - dists.min())
pltr = pv.Plotter(window_size=[512, 512], off_screen=False)
pltr.background_color = "#363940"
pltr.set_focus(mesh.center)
pltr.set_position((17, 14, 12))
pltr.camera.zoom(1)
mesh["dist"] = dists
pltr.add_mesh(
    mesh,
    smooth_shading=True,
    specular=0.9,
    color="yellow
)
pltr.show()

I want to clip it to the ball of radius 3. That is to say, I want to discard all points whose distance to the origin is superior to 3. How could I do? The boundary where the mesh is cut must be smooth.

PS: it seems PyVista has changed since the last time I used it, I'm not sure everything in my code is needed.

Upvotes: 0

Views: 221

Answers (1)

Stéphane Laurent
Stéphane Laurent

Reputation: 84659

Ok, I get it, using spherical coordinates:

import numpy as np
import pyvista as pv

def f(ρ, θ, ϕ, A, B):
    x = ρ * np.cos(θ) * np.sin(ϕ)
    y = ρ * np.sin(θ) * np.sin(ϕ)
    z = ρ * np.cos(ϕ)    
    return (
        z**4 * B**2
        + 4 * x * y**2 * A * B**2
        + x * z**2 * A * B**2
        - 2 * z**4 * A
        - 4 * x * y**2 * B**2
        - x * z**2 * B**2
        + 3 * z**2 * A * B**2
        - 2 * z**4
        - x * A * B**2
        - 2 * z**2 * A
        + x * B**2
        + A * B**2
        + 2 * z**2
        - B**2
    )

# generate data grid for computing the values
Rho, Theta, Phi = np.mgrid[0:np.sqrt(3):140j, 0:(2*np.pi):140j, 0:np.pi:140j]
A = np.cos(2.5*np.pi/4)
B = np.sin(2.5*np.pi/4)
values = f(Rho, Theta, Phi, A, B)

# create a structured grid
mesh = pv.StructuredGrid(Rho, Theta, Phi)
mesh.point_data['values'] = values.ravel(order='F')  # also the active scalars

def sph2cart(sph):
    ρ = sph[0]
    θ = sph[1]
    ϕ = sph[2]
    return [
        ρ * np.cos(θ) * np.sin(ϕ),
        ρ * np.sin(θ) * np.sin(ϕ),
        ρ * np.cos(ϕ)
    ]
mesh.points = np.apply_along_axis(sph2cart, 1, mesh.points)

# compute one isosurface
isos = mesh.contour(isosurfaces=1, rng=[0,0])

# plot it interactively 
isos.plot(smooth_shading=True, show_scalar_bar=False)

But this is applicable only for an isosurface. I'm still wondering how to do for a parametric mesh.

Upvotes: 1

Related Questions