svVascularize

svv.domain.Domain

Core Class v1.0.0

The Domain class defines the region in space that will be recognized by svVascularize when generating vascular networks. The class abstracts the physical representation of the space to allow for efficient interrogation and data manipulation.

Domains can be created from point clouds, PyVista objects, or by combining existing domains using set operations. The class provides implicit function evaluation, boundary extraction, and interior mesh generation capabilities.

Quick Example

import numpy as np
import pyvista as pv
from svv.domain.domain import Domain

# Create domain from PyVista mesh
mesh = pv.read('vessel_geometry.stl')
domain = Domain(mesh)

# Build implicit function
domain.create()
domain.solve()
domain.build(resolution=25)

# Get boundary and interior mesh
boundary, grid = domain.get_boundary(resolution=30)
interior_mesh = domain.get_interior()

# Sample points from domain
interior_points, _ = domain.get_interior_points(n=1000)
boundary_points = domain.get_boundary_points(n=500)

Constructor

Domain(*args, **kwargs)

Parameters

Parameter Type Description
*args variable Can be:
  • Single numpy array of points (n×d)
  • Two numpy arrays: points (n×d) and normals (n×d)
  • PyVista object (PolyData, UnstructuredGrid, etc.)
**kwargs dict Additional keyword arguments passed to data processing methods

Attributes

Attribute Type Description
points np.ndarray Point cloud defining the domain (n×d array)
normals np.ndarray Normal vectors at each point (n×d array)
boundary pv.PolyData Boundary mesh of the domain
mesh pv.UnstructuredGrid Interior tetrahedral/triangular mesh
patches list List of Patch objects for implicit function
functions list List of implicit functions
n int Number of points
d int Spatial dimension (2 or 3)
characteristic_length float Characteristic length scale of the domain
area float Surface area (3D) or perimeter (2D)
volume float Volume (3D) or area (2D)
convexity float Ratio of domain volume to convex hull volume

Methods

Core Methods

set_data(*args, **kwargs)

Set the data for the domain from point-wise data or a PyVista object.

Parameters

  • *args: Points array, points and normals arrays, or PyVista object
  • **kwargs: Additional processing parameters
create(min_patch_size=10, max_patch_size=20, overlap=0.2, feature_angle=30)

Create patches for the domain. This is the first step in building the implicit function.

Parameters

  • min_patch_size (int, default=10): Minimum number of points required to form a patch.
  • max_patch_size (int, default=20): Target maximum number of points per patch (nearest-neighbor window).
  • overlap (float, default=0.2): Fraction [0, 1] controlling allowed overlap of points between patches.
  • feature_angle (float, default=30): Maximum angle in degrees between point-wise normal vectors for inclusion in the same patch (used only when normals are provided).
solve(method="L-BFGS-B", precision=9)

Solve the individual patch interpolation problems prior to blending.

Parameters

  • method (str, default="L-BFGS-B"): Optimization method for SciPy's minimize (see Patch.solve).
  • precision (int, default=9): Decimal places to round the solved constants.
build(resolution=25, skip_boundary=False)

Build the implicit function describing the domain and optionally extract boundary/mesh artifacts.

Parameters

  • resolution (int, default=25): Grid resolution for boundary extraction.
  • skip_boundary (bool, default=False): If true, only assemble fast-evaluation structures and skip boundary and interior mesh generation.

Evaluation Methods

__call__(points, **kwargs)

Evaluate the implicit function at given points.

Parameters

  • points (np.ndarray): Points to evaluate (m×d array)
  • k (int, default=1): Number of nearest patches to consider
  • normalize (bool, default=True): Normalize output values
  • tolerance (float): Tolerance for patch selection

Returns

  • values (np.ndarray): Implicit function values (m×1 array)
within(points, level=0, **kwargs)

Determine if points are within the domain.

Returns

  • inside (np.ndarray): Boolean array indicating if points are inside

Mesh Generation

get_boundary(resolution, **kwargs)

Extract the boundary surface of the domain.

Parameters

  • resolution (int): Grid resolution for marching cubes/squares
  • get_largest (bool, default=True): Reserved (current implementation always returns the largest connected component)

Returns

  • boundary (pv.PolyData): Boundary mesh
  • grid (pv.ImageData): Volumetric grid used for extraction
get_interior(verbose=False, **kwargs)

Generate tetrahedral (3D) or triangular (2D) mesh of the interior.

Parameters

  • verbose (bool): Print mesh generation progress
  • **kwargs: Parameters passed to TetGen/Triangle

Returns

  • mesh (pv.UnstructuredGrid): Interior mesh

Point Sampling

get_interior_points(n, tree=None, volume_threshold=None, threshold=None, method=None, implicit_range=(-1.0, 0.0), **kwargs)

Sample random points from the domain interior.

Parameters

  • n (int): Number of points to sample
  • tree (KDTree, optional): Existing tree structure for collision detection
  • volume_threshold (float, optional): Maximum distance from tree
  • threshold (float, optional): Minimum distance from tree
  • method (str, optional): Sampling method ('implicit_only', 'preallocate', or default)
  • implicit_range (tuple): Range of implicit function values to accept

Returns

  • points (np.ndarray): Sampled points (n×d array)
  • cells (np.ndarray): Cell indices for each point
get_boundary_points(n, method=None, **kwargs)

Sample random points from the domain boundary.

Parameters

  • n (int): Number of points to sample
  • method (str, optional): Sampling method

Returns

  • points (np.ndarray): Sampled boundary points (n×d array)

Utility Methods

set_random_seed(seed)

Set the random seed for reproducible sampling.

set_random_generator()

Initialize the random number generator with the set seed.

Operators

domain1 + domain2

Union of two domains (set union operation).

# Create union of two domains
combined_domain = domain1 + domain2
domain1 - domain2

Difference of two domains (set difference operation).

# Subtract domain2 from domain1
difference_domain = domain1 - domain2

Examples

Creating a Domain from STL File

import pyvista as pv
from svv.domain.domain import Domain

# Load STL file
mesh = pv.read('vessel.stl')

# Create domain
domain = Domain(mesh)

# Build implicit function with custom resolution
domain.create()
domain.solve()
domain.build(resolution=30)

# Get mesh statistics
print(f"Domain volume: {domain.volume:.2f}")
print(f"Domain surface area: {domain.area:.2f}")
print(f"Convexity: {domain.convexity:.3f}")

Customized Patch Allocation

import pyvista as pv
from svv.domain.domain import Domain

mesh = pv.read('vessel.stl')
domain = Domain(mesh)

# Tune patch allocation for dense data with noisy normals
domain.create(
    min_patch_size=8,      # allow smaller patches in sparse areas
    max_patch_size=32,     # include more neighbors in dense regions
    overlap=0.35,          # increase overlap for smoother blending
    feature_angle=45       # be more tolerant to normal variation/noise
)
domain.solve()
domain.build(resolution=30)
When to adjust defaults:
  • Dense point clouds: increase max_patch_size to capture local context and reduce edge artifacts.
  • Sparse/heterogeneous sampling: decrease min_patch_size so isolated regions still form patches.
  • Noisy or inconsistent normals: raise feature_angle to admit neighbors with larger normal differences; lower it to preserve sharp features.
  • Smoother blends vs. speed: higher overlap yields smoother transitions across patches but increases compute.

Boolean Operations on Domains

# Create two domains
sphere1 = pv.Sphere(center=(0, 0, 0), radius=1.0)
domain1 = Domain(sphere1)
domain1.create()
domain1.solve()
domain1.build()

sphere2 = pv.Sphere(center=(0.5, 0, 0), radius=0.8)
domain2 = Domain(sphere2)
domain2.create()
domain2.solve()
domain2.build()

# Create union (merged domains)
union_domain = domain1 + domain2

# Create difference (domain1 with hole from domain2)
diff_domain = domain1 - domain2

# Evaluate at test points
test_points = np.random.randn(100, 3)
union_values = union_domain(test_points)
diff_values = diff_domain(test_points)

Sampling Points with Collision Avoidance

from scipy.spatial import cKDTree

# Create domain
domain = Domain(mesh)
domain.create()
domain.solve()
domain.build()

# Create existing tree structure (e.g., from vessel centerlines)
existing_points = np.array([[0, 0, 0], [1, 0, 0], [0, 1, 0]])
tree = cKDTree(existing_points)

# Sample points avoiding existing structure
new_points, cells = domain.get_interior_points(
    n=1000,
    tree=tree,
    threshold=0.5,  # Minimum distance from existing points
    volume_threshold=5.0  # Maximum distance from existing points
)

# Sample boundary points
boundary_points = domain.get_boundary_points(n=500)

Notes

Implicit Function: The domain uses radial basis function (RBF) interpolation to create a smooth implicit function representation. Points inside the domain have negative values, points outside have positive values, and the boundary is at the zero level set.

Memory Usage: For large point clouds or high-resolution meshes, the implicit function construction can be memory-intensive. Consider using lower resolution for initial testing.

Performance: The domain uses optimized Cython routines for point sampling (c_sample.pyx) and allocation (c_allocate.pyx). These provide significant speedups for large-scale operations.

Best Practice: Always call methods in order: create() → solve() → build() to properly construct the implicit function before using evaluation or sampling methods.

See Also