Source code for homcloud.interface.bitmap_phtrees

import numpy as np

import homcloud.pict.tree as pict_tree
import homcloud.paraview_interface as pv_interface
import homcloud.plotly_3d as p3d
from homcloud.spatial_searcher import SpatialSearcher
from . import bitmap as bm


[docs] class BitmapPHTrees(object): """ This class represents PH trees computed from a bitmap filtration. You can create an instance of this class by :meth:`PDList.bitmap_phtrees`. Please see `Obayashi (2018) <https://doi.org/10.1137/17M1159439>`_ if you want to know more about optimal volumes, and see `Obayashi (2023) <https://doi.org/10.1007/s41468-023-00119-8>`_ if you want to know more about stable volumes. Attributes: nodes (list of :class:`BitmapPHTrees.Node`): All nodes essential_nodes (list of :class:`BitmapPHTrees.Node`): All essential nodes nonessential_nodes (list of :class:`BitmapPHTrees.Node`): All nonessential nodes roots (list of :class:`BitmapPHTrees.Node`): All root nodes """
[docs] @staticmethod def for_bitmap_levelset(array, mode="sublevel", save_to=None): """Computes a PH trees from 0-th and (n-1)-st persistent homology from levelset filtration of the bitmap. Args: array (numpy.ndarray): The bitmap data. mode (string): "superlevel" or "sublevel". save_to (string or None): The filepath which the PH trees is stored in. Returns: :class:`PDList`: The 0th and (n-1)st PDs with PH-trees. """ from .pd import PDList is_superlevel = mode == "superlevel" lower, upper = pict_tree.construct_mergetrees(array, is_superlevel) f = PDList.open_pdgm_file(save_to) pict_tree.save_pdgm(f, array.ndim, is_superlevel, lower, upper) f.seek(0) return PDList(f, PDList.FileType.PDGM)
def __init__(self, treedict, is_superlevel): self.treedict = treedict self.is_superlevel = is_superlevel self.nodes = [self.node(id) for id in self.treedict["nodes"] if not self.is_trivial_id(id)] self.essential_nodes = [node for node in self.nodes if node.essential()] self.nonessential_nodes = [node for node in self.nodes if not node.essential()] self.spatial_searcher = self.build_spatial_searcher() self.roots = [node for node in self.nodes if node.isroot()] def node(self, id): return BitmapPHTrees.Node(self, self.treedict["nodes"][id]) def is_trivial_id(self, id): node = self.treedict["nodes"][id] return node["birth-time"] == node["death-time"] def build_spatial_searcher(self): births = np.array([node.birth_time() for node in self.nonessential_nodes]) deaths = np.array([node.death_time() for node in self.nonessential_nodes]) return SpatialSearcher(self.nonessential_nodes, births, deaths)
[docs] def degree(self): """ Returns: int: The degree of the persistent homology. """ return self.treedict["degree"]
[docs] def nearest_pair_node(self, x, y): """ Searches a tree node corresponding the birth-death pair nearest to (x, y). Args: x (float): The x-coordinate. y (float): The y-coordinate. Returns: :class:`BitmapPHTrees.Node`: The node. """ return self.spatial_searcher.nearest_pair(x, y)
[docs] def pair_nodes_in_rectangle(self, xmin, xmax, ymin, ymax): """ Searches tree nodes corresponding the birth-death pairs which the given rectangle contains. Args: xmin (float): The left side of the rectangle. xmax (float): The right side of the rectangle. ymin (float): The bottom side of the rectangle. ymax (float): The top side of the rectangle. Returns: list of :class:`BitmapPHTrees.Node`: The nodes. """ return self.spatial_searcher.in_rectangle(xmin, xmax, ymin, ymax)
def sign(self): return -1 if self.is_superlevel else 1
[docs] class Volume(object): """ This class is the superclass of :class:`Node` and :class:`StableVolume`. This class defines common methods for these two classes. """
[docs] def to_paraview_node(self, gui_name=None): """ Construct a :class:`homcloud.paraview_interface.PipelineNode` object to visulize :meth:`volume`. You can show the volume by :meth:`homcloud.paraview_interface.show`. You can also adjust the visual by the methods of :class:`homcloud.paraview_interface.PipelineNode`. Args: gui_name (string or None): The name shown in Pipeline Browser in paraview's GUI. Returns: homcloud.paraview_interface.PipelineNode: A PipelineNode object. """ return BitmapPHTrees.to_paraview_node_from_nodes([self], gui_name)
to_pvnode = to_paraview_node
[docs] def to_plotly3d_trace(self, color=None, name=""): """ Constructs a plotly's trace object to visualize the volume Args: color (string or None): The name of the color name (string): The name of the object Returns: plotly.graph_objects.Mesh3d: Plotly's trace object """ return p3d.Voxels(self.volume(), color, name)
to_plotly3d = to_plotly3d_trace
[docs] def to_pyvista_mesh(self): """ Constructs a PyVista's mesh object to visualize an optimal/stable volume Returns: pyvista.PolyData: PyVista's mesh object """ import homcloud.pyvistahelper as pvhelper return pvhelper.SparseVoxels(self.volume())
[docs] class Node(Volume): """ This class represents a tree node of PH trees for a bitmap filtration. You can draw the optimal volume corresponding to the node on an image by :meth:`draw_volumes_on_2d_image`. """ def __init__(self, mt, dic): self.mt = mt self.dic = dic
[docs] def birth_time(self): """ Returns: float: The birth time. """ return self.dic["birth-time"]
[docs] def death_time(self): """ Returns: float: The death time. May be infinity. """ if not self.essential(): return self.dic["death-time"] return self.mt.sign() * np.inf
[docs] def lifetime(self): """ Returns: float: The lifetime of the pair. """ return self.death_time() - self.birth_time()
def __iter__(self): return (self.birth_time(), self.death_time()).__iter__()
[docs] def birth_position(self): """Returns the birth pixel. Returns: list of int: The coordinate of the birth pixel. """ return self.dic["birth-pixel"]
birth_pixel = birth_position
[docs] def death_position(self): """Returns the death pixel. Returns: list of int: The coordinate of the death pixel. """ return self.dic["death-pixel"]
death_pixel = death_position
[docs] def volume(self): """ Returns the optimal volume. In fact, for degree 0 node, this method returns optimal 0-cohomological volume (this means maximal connected conponents in the filtration) and for degree (n-1), this method returns optimal (n-1)-homological volume. This fact is quite mathematical problem, so if you feel that it is too difficult to understand, you can ignore the fact. You can understand the volume information without the understanding of such mathematical background. Returns: list of list of int: The coordinates of all pixels in the optimal volume. """ return self.dic["volume"]
[docs] def stable_volume(self, epsilon): """ Returns the stable volume of the optimal volume with parameter epsilon. Args: epsilon (float): Duration noise strength Returns: :class:`BitmapPHTrees.StableVolume`: The stable volume """ return BitmapPHTrees.StableVolume(self, epsilon)
[docs] def essential(self): """ Returns: bool: True if the death time is infinity. """ return self.dic["death-time"] is None
[docs] def parent(self): """ Returns: BitmapPHTrees.Node: The parent node of the node in the PH tree. """ parent_id = self.dic["parent"] return BitmapPHTrees.Node(self.mt, self.mt.treedict["nodes"][parent_id])
[docs] def isroot(self): """ Returns: bool: True if the node is a root node. """ return self.dic["parent"] is None
[docs] def children(self): """ Returns: list of :class:`BitmapPHTrees.Node`: All children nodes. """ return [self.mt.node(child_id) for child_id in self.dic["children"] if not self.mt.is_trivial_id(child_id)]
def __eq__(self, other): return (self.mt == other.mt) and (self.dic["id"] == other.dic["id"]) def __repr__(self): return "BitmapPHTrees.Node({}, {})".format(self.birth_time(), self.death_time())
[docs] class StableVolume(Volume): """ This class represents a stable volume in :class:`BitmapPHTrees`. """ def __init__(self, node, epsilon): self.node = node self._volume = set(tuple(p) for p in node.dic["volume"]) for child_id in node.dic["children"]: child_dic = node.mt.treedict["nodes"][child_id] if self.is_unstable_child(child_dic, epsilon): self._volume.difference_update(tuple(p) for p in child_dic["volume"]) def is_unstable_child(self, child_dic, epsilon): is_sublevel = not self.node.mt.is_superlevel is_lower = self.node.mt.degree() == 0 if is_lower and is_sublevel: return child_dic["death-time"] > self.node.death_time() - epsilon elif is_lower and not is_sublevel: return child_dic["death-time"] < self.node.death_time() + epsilon elif not is_lower and is_sublevel: return child_dic["birth-time"] < self.node.birth_time() + epsilon elif not is_lower and not is_sublevel: return child_dic["birth-time"] > self.node.birth_time() - epsilon
[docs] def volume(self): """ Returns the pixels of the stable volume. Returns: list[tuple[int,...]]: The coordinates of all pixels. """ return list(self._volume)
[docs] @staticmethod def to_paraview_node_from_nodes(nodes, gui_name=None): """ Construct a :class:`homcloud.paraview_interface.PipelineNode` object to visulize node volumes of the nodes. Args: nodes (list of :class:`BitmapPHTrees.Node`): The list of nodes to be visualized. gui_name (string or None): The name shown in Pipeline Browser in paraview's GUI. Returns: homcloud.paraview_interface.VTK: Paraview Pipeline node object. """ all_volume = [] for node in nodes: all_volume.extend(node.volume()) min, max = bm.volume_range(all_volume) bitmap = bm.volume_to_bitmap(all_volume, (min, max)) return pv_interface.VoxelData(bitmap, gui_name, min).threshold("value", (0.5, 1.0))
to_pvnode_from_nodes = to_paraview_node_from_nodes