API Reference

Automatically generated documentation for core functions in pyaof library.

Core Functions

Compute Average Outward Flux

`pyaof.compute_aof(sdf_np, hx=1.0, hy=1.0, hz=1.0, tile_size=128, halo=2)`

Compute Average Outward Flux (AOF) from a Signed Distance Field.

This function computes the Average Outward Flux (AOF) volume from a Signed Distance Field (SDF) and applies local non-maximum suppression to enhance ridge-like structures and suppress weak responses.

Parameters:

Name	Type	Description	Default
`sdf_np`	`ndarray`	Input Signed Distance Field (SDF) volume as a 2D/3D NumPy array.	required
`hx`	`float`	Grid spacing along the x-axis. Defaults to `1.0`.	`1.0`
`hy`	`float`	Grid spacing along the y-axis. Defaults to `1.0`.	`1.0`
`hz`	`float`	Grid spacing along the z-axis. Defaults to `1.0`.	`1.0`
`tile_size`	`int`	Tile size used during tiled AOF computation. Defaults to `128`. Larger tile sizes may improve performance on larger GPUs but increase memory usage.	`128`
`halo`	`int`	Halo padding size used for neighborhood access during tiled computation. Defaults to `2`.	`2`

Returns:

Type	Description
`ndarray`	np.ndarray: A thresholded 3D AOF volume after non-maximum suppression.

Example

aof_vol = compute_aof(sdf)

Typical Use Cases

Medial axis approximation
Skeleton extraction

`pyaof.quick_aof(sdf_np, hx=1.0, hy=1.0, hz=1.0, tile_size=128, halo=2)`

Computes the Average Outward Flux (AOF) map for a 2D or 3D Signed Distance volume.

Parameters:

Name	Type	Description	Default
`sdf_np`	`ndarray`	A 3D (Z, Y, X) or 2D (Y, X) NumPy array representing the Signed Distance Field or intensity volume.	required

Returns: aof (ndarray): A NumPy array of the same shape as 'sdf_np'.

Compute signed distance function

`pyaof.compute_sdf(mesh, resolution, padding=0.1, normalize=False, smoothing=True)`

Compute a Signed Distance Field (SDF) from a VTK mesh.

This function voxelizes the input mesh, builds an occupancy grid, and computes a signed distance field (SDF) using ITK signed distance transform.

Parameters:

Name	Type	Description	Default
`mesh`	`vtkPolyData`	Input VTK mesh.	required
`resolution`	`int`	Resolution of the voxel grid along the longest axis.	required
`padding`	`float`	Extra padding around mesh bounds to avoid clipping artifacts. Defaults to 0.1.	`0.1`
`normalize`	`bool`	If True, scales mesh into a unit cube before voxelization. Defaults to False.	`False`
`smoothing`	`bool`	If True, smoothes the sdf to avoid staircasing problems at the boundary. Defaults to True.	`True`

Returns:

Type Description

tuple[ndarray, vtkPolyData, ndarray, tuple[float, float, float]]

tuple[np.ndarray, vtk.vtkPolyData, np.ndarray, tuple[float, float, float]]: A tuple containing:

sdf (np.ndarray): Signed distance field volume of shape (N, N, N).
normalized_mesh (vtk.vtkPolyData): Mesh used for voxelization (normalized if enabled).
coords (np.ndarray): 3D grid coordinates corresponding to voxel centers.
spacing (tuple[float, float, float]): Grid spacing in (x, y, z) directions.

Notes

SDF sign convention: negative inside, positive outside.
Requires VTK for mesh processing.
Distance computation uses Euclidean distance transform.

Example

sdf, mesh_n, coords, spacing = compute_sdf( ... mesh, ... resolution=1024, ... normalize=True ... )

Utilities

Meshing Utilities

`pyaof.voxelize_vtk(mesh, spacing=(1.0, 1.0, 1.0), origin=(-0.5, -0.5, -0.5), resolution=256, padding=0.1, inside_value=1, outside_value=0, normalize=False)`

Voxelize a VTK mesh into a 3D binary volume.

This function converts a triangular VTK mesh into a structured 3D voxel grid using VTK's image stencil pipeline. It performs several preprocessing steps (cleaning, triangulation, hole filling, and normal correction) before rasterizing the mesh into a volumetric grid.

Parameters:

Name	Type	Description	Default
`mesh`	`vtkPolyData`	Input triangular surface mesh.	required
`spacing`	`tuple[float, float, float]`	Voxel spacing in (x, y, z). Defaults to (1.0, 1.0, 1.0).	`(1.0, 1.0, 1.0)`
`origin`	`tuple[float, float, float]`	Origin of the voxel grid. Defaults to (-0.5, -0.5, -0.5).	`(-0.5, -0.5, -0.5)`
`resolution`	`int`	Number of voxels along each axis (NxNxN grid). Defaults to 256.	`256`
`padding`	`float`	Extra padding applied during normalization to avoid clipping. Defaults to 0.1.	`0.1`
`inside_value`	`int`	Value assigned to voxels inside the mesh. Defaults to 1.	`1`
`outside_value`	`int`	Value assigned to voxels outside the mesh. Defaults to 0.	`0`
`normalize`	`bool`	If True, scales mesh into a unit cube before voxelization and overrides spacing/origin accordingly. Defaults to False.	`False`

Returns:

Name	Type	Description
`tuple`	`tuple[ndarray, tuple[float, float, float], vtkPolyData, ndarray]`	A tuple containing: np.ndarray: 3D voxel grid of shape (Nz, Ny, Nx) containing binarized object. tuple[float, float, float]: Effective voxel spacing used for the grid. vtk.vtkPolyData: Processed (cleaned + normalized) mesh used for voxelization. np.ndarray: Vertex coordinates of the processed mesh in reordered (z, y, x) format.

Example

vol, spacing, mesh_n, coords = voxelize_vtk( ... mesh, ... resolution=1024, ... normalize=True ... )

`pyaof.normalize_mesh_to_unit(mesh, padding=0.01)`

Normalize a VTK mesh into a unit cube centered at the origin. A small padding is added to avoid boundary clipping.

Parameters:

Name	Type	Description	Default
`mesh`	`vtkPolyData`	Input triangular mesh to normalize.	required
`padding`	`float`	Extra margin added to the mesh bounding radius before scaling. Prevents tight clipping at the unit cube boundary. Defaults to 0.01.	`0.01`

Returns:

Type	Description
`vtkPolyData`	vtk.vtkPolyData: Normalized mesh centered at the origin and scaled to fit inside a unit cube.

Example

normalized = normalize_mesh_to_unit(mesh)

Point Sampling

`pyaof.sample_sphere_points(n_points, *, dimension=3, iterations=500)`

Generate approximately uniform points on a unit sphere. This function generates points on the surface of a unit sphere using an iterative electrostatic repulsion method. Points are initialized randomly and repeatedly pushed apart using inverse-square repulsive forces to improve uniformity.

The first point is fixed at [1, 0, 0, ...] to stabilize the optimization and remove rotational ambiguity.

Parameters:

Name	Type	Description	Default
`n_points`	`int`	Number of points to generate on the sphere surface.	required
`dimension`	`int`	Dimensionality of the embedding space. Defaults to `3`. `dimension=3` generates points on a standard 3D sphere. Higher values generate points on hyperspheres.	`3`
`iterations`	`int`	Number of relaxation iterations used to improve point distribution. Defaults to `500`.	`500`

Returns:

Type	Description
`ndarray`	np.ndarray: Array of shape `(n_points, dimension)` containing normalized direction vectors on the unit sphere surface.

Example

dirs = sample_sphere_points( ... n_points=128, ... dimension=3, ... iterations=100 ... ) dirs.shape (128, 3)

Source code in pyaof/aof.py

def sample_sphere_points(n_points:int, *, dimension=3, iterations=500)->np.ndarray:
    """Generate approximately uniform points on a unit sphere.
    This function generates points on the surface of a unit sphere using
    an iterative electrostatic repulsion method. Points are initialized
    randomly and repeatedly pushed apart using inverse-square repulsive
    forces to improve uniformity.

    The first point is fixed at ``[1, 0, 0, ...]`` to stabilize the
    optimization and remove rotational ambiguity.

    Args:
        n_points (int):
            Number of points to generate on the sphere surface.

        dimension (int, optional):
            Dimensionality of the embedding space. Defaults to ``3``.

            - ``dimension=3`` generates points on a standard 3D sphere.
            - Higher values generate points on hyperspheres.

        iterations (int, optional):
            Number of relaxation iterations used to improve point
            distribution. Defaults to ``500``.

    Returns:
        np.ndarray:
            Array of shape ``(n_points, dimension)`` containing normalized
            direction vectors on the unit sphere surface.

    Example:
        >>> dirs = sample_sphere_points(
        ...     n_points=128,
        ...     dimension=3,
        ...     iterations=100
        ... )
        >>> dirs.shape
        (128, 3)

    """
    dirs = np.random.rand(n_points, dimension)
    # First point is fixed at [1, 0, 0...], others are random
    dirs[0, :] = 0
    dirs[0, 0] = 1.0
    dirs /= np.linalg.norm(dirs, axis=-1, keepdims=True)

    for _ in range(iterations):
        # 2. Calculate pairwise difference vectors: (N, N, Dim)
        # points[:, np.newaxis, :] is (N, 1, Dim)
        # points[np.newaxis, :, :] is (1, N, Dim)
        diff = dirs[:, np.newaxis, :] - dirs[np.newaxis, :, :]

        r2 = np.sum(diff**2, axis=-1)
        np.fill_diagonal(r2, np.inf)
        forces_ij = diff / r2[:, :, np.newaxis]
        total_forces = np.sum(forces_ij, axis=1)

        # Points[0] is fixed, so we only update points[1:] and Re-normalize
        dirs[1:] += total_forces[1:]
        dirs[1:] /= np.linalg.norm(dirs[1:], axis=-1, keepdims=True)

    return dirs

Visualize 360 rotation

`pyaof.util.render_360_rotation(target=np.array([0.0, 0.0, 0.0]), radius=5.0, elevation=2.0, num_frames=60, output_folder='renders', up_axis='Z')`

Renders a 360-degree camera rotation around a target point in Polyscope. Safely exits if Polyscope is not installed or initialized.

Parameters:

Name	Type	Description	Default
`target`	`ndarray`	The 3D coordinate [X, Y, Z] the camera will focus on. Defaults to np.array([0.0, 0.0, 0.0]).	`array([0.0, 0.0, 0.0])`
`radius`	`float`	The horizontal distance from the camera to the target point. Defaults to 5.0.	`5.0`
`elevation`	`float`	The vertical height offset of the camera relative to the target's height. Defaults to 2.0.	`2.0`
`num_frames`	`int`	The total number of subdivisions along the circular camera path. Defaults to 60.	`60`
`up_axis`	`str`	The global vertical vector convention of the dataset. Must be either 'X','Y' or 'Z' (case-insensitive). Defaults to "Z".	`'Z'`

Raises:

Type	Description
`ValueError`	If `up_axis` is passed as anything other than 'Y' or 'Z'.

Returns:

Type	Description
	None

Examples:

>>> import numpy as np
>>> import polyscope as ps
>>> ps.init()
>>> # ... register your point clouds or meshes here ...
>>> render_360_rotation(target=np.array([0, 1, 0]), radius=3.5, up_axis="Y")