Namespace hg¶
Warning
The cpp documentation is in progress, so here is a raw dump of the whole hg
namespace: enjoy!
-
enum class
hg
::
accumulators
¶ Values:
-
enumerator
first
¶
-
enumerator
last
¶
-
enumerator
mean
¶
-
enumerator
min
¶
-
enumerator
max
¶
-
enumerator
counter
¶
-
enumerator
sum
¶
-
enumerator
prod
¶
-
enumerator
argmin
¶
-
enumerator
argmax
¶
-
enumerator
-
enum class
hg
::
weight_functions
¶ Predefined edge-weighting functions (see weight_graph function)
Values:
-
enumerator
mean
¶
-
enumerator
min
¶
-
enumerator
max
¶
-
enumerator
L0
¶
-
enumerator
L1
¶
-
enumerator
L2
¶
-
enumerator
L_infinity
¶
-
enumerator
L2_squared
¶
-
enumerator
source
¶
-
enumerator
target
¶
-
enumerator
-
enum class
hg
::
optimal_cut_measure
¶ Values:
-
enumerator
BCE
¶
-
enumerator
DHamming
¶
-
enumerator
DCovering
¶
-
enumerator
-
enum class
hg
::
partition_measure
¶ Values:
-
enumerator
BCE
¶
-
enumerator
DHamming
¶
-
enumerator
DCovering
¶
-
enumerator
-
enum
hg
::
tos_padding
¶ Padding mode for the function component_tree_tree_of_shapes.
Values:
-
enumerator
none
¶
-
enumerator
mean
¶
-
enumerator
mean
-
enumerator
mean
-
enumerator
zero
¶
-
enumerator
-
enum class
hg
::
leaves_it
¶ Enum used in tree node iterator (leaves_to_root_iterator and root_to_leaves_iterator) to include or exclude leaves from the iterator.
Values:
-
enumerator
include
¶
-
enumerator
exclude
¶
-
enumerator
-
enum class
hg
::
root_it
¶ Enum used in tree node iterator (leaves_to_root_iterator and root_to_leaves_iterator) to include or exclude the root from the iterator.
Values:
-
enumerator
include
¶
-
enumerator
exclude
¶
-
enumerator
-
using
hg
::
contour_2d
= contour_2d_internal::contour_2d<point_2d_f>¶
-
using
hg
::
polyline_contour_2d
= contour_2d_internal::polyline_contour_2d<point_2d_f>¶
-
using
hg
::
contour_segment_2d
= contour_2d_internal::contour_segment_2d<point_2d_f>¶
-
template<int
dim
>
usinghg
::
embedding_grid
= typename embedding_internal::embedding_grid<dim, index_t>¶
-
template<typename
value_type
>
usinghg
::
fibonacci_heap
= fibonacci_heap_internal::fibonacci_heap<value_type>¶
-
using
hg
::
lca_sparse_table_block
= lca_internal::lca_rmq<tree, range_minimum_query_internal::rmq_sparse_table_block<index_t>>¶
-
using
hg
::
lca_sparse_table
= lca_internal::lca_rmq<tree, range_minimum_query_internal::rmq_sparse_table<index_t>>¶
-
using
hg
::
lca_fast
= lca_sparse_table_block¶
-
template<typename
value_t
, unsigned intdim
>
usinghg
::
point
= xt::xtensor_fixed<value_t, xt::xshape<dim>>¶
-
template<typename
embedding_t
>
usinghg
::
regular_graph
= regular_graph_internal::regular_graph<embedding_t>¶
-
using
hg
::
regular_grid_graph_1d
= regular_graph<hg::embedding_grid_1d>¶
-
using
hg
::
regular_grid_graph_2d
= regular_graph<hg::embedding_grid_2d>¶
-
using
hg
::
regular_grid_graph_3d
= regular_graph<hg::embedding_grid_3d>¶
-
using
hg
::
regular_grid_graph_4d
= regular_graph<hg::embedding_grid_4d>¶
-
template<typename
embedding_t
>
usinghg
::
regular_graph_out_edge_iterator
= typename regular_graph_internal::regular_graph<embedding_t>::out_edge_iterator¶
-
template<typename
embedding_t
>
usinghg
::
regular_graph_adjacent_vertex_iterator
= typename regular_graph_internal::regular_graph<embedding_t>::adjacency_iterator¶
-
using
hg
::
tree
= tree_internal::tree¶
-
using
hg
::
vecS
= undirected_graph_internal::vecS¶
-
using
hg
::
hash_setS
= undirected_graph_internal::hash_setS¶
-
template<typename
storage_type
= vecS>
usinghg
::
undirected_graph
= undirected_graph_internal::undirected_graph<storage_type>¶
-
using
hg
::
ugraph
= undirected_graph_internal::undirected_graph<>¶
-
using
hg
::
union_find
= union_find_internal::union_find<>¶
-
using
hg
::
index_t
= int64_t¶ Preferred type to represent an index.
-
using
hg
::
size_t
= std::size_t¶ Preferred type to represent a size.
-
const index_t
hg
::
invalid_index
= -1¶ Constant used to represent an invalid index (eg.
not initialized)
-
template<typename
T
, typenameaccumulator_t
, typenameoutput_t
= typename T::value_type>
autohg
::
accumulate_at
(const array_1d<index_t> &indices, const xt::xexpression<T> &xweights, const accumulator_t &accumulator)¶ Accumulate the given weights located at given indices.
Let :math:
M = max(indices)
. For all :math:i \in \{0, \ldots, M\}
.. math::
result[i] = accumulator(\{weights[j, :] \mid indices[j] = i \})
- Template Parameters:
T –
accumulator_t –
output_t –
- Parameters:
indices – a 1d array of indices (entry equals to :math:
-1
are ignored)xweights – a nd-array of shape :math:
(s_1, \ldots, s_n)
such that :math:s_1=indices.size()
accumulator –
- Returns:
a nd-array of size :math:
(M, s_2, \ldots, s_n)
-
template<typename
graph_t
, typenameT
, typenameaccumulator_t
, typenameoutput_t
= typename T::value_type>
autohg
::
accumulate_graph_edges
(const graph_t &graph, const xt::xexpression<T> &xedge_weights, const accumulator_t &accumulator)¶
-
template<typename
graph_t
, typenameT
, typenameaccumulator_t
, typenameoutput_t
= typename T::value_type>
autohg
::
accumulate_graph_vertices
(const graph_t &graph, const xt::xexpression<T> &xvertex_weights, const accumulator_t &accumulator)¶
-
template<typename
tree_t
, typenameT
, typenameaccumulator_t
, typenameoutput_t
= typename T::value_type>
autohg
::
accumulate_parallel
(const tree_t &tree, const xt::xexpression<T> &xinput, const accumulator_t &accumulator)¶
-
template<typename
tree_t
, typenameT
, typenameaccumulator_t
, typenameoutput_t
= typename T::value_type>
autohg
::
accumulate_sequential
(const tree_t &tree, const xt::xexpression<T> &xvertex_data, const accumulator_t &accumulator)¶
-
template<typename
tree_t
, typenameT1
, typenameT2
, typenameaccumulator_t
, typenamecombination_fun_t
, typenameoutput_t
= typename T1::value_type>
autohg
::
accumulate_and_combine_sequential
(const tree_t &tree, const xt::xexpression<T1> &xinput, const xt::xexpression<T2> &xvertex_data, const accumulator_t &accumulator, const combination_fun_t &combine)¶
-
template<typename
tree_t
, typenameT1
>
autohg
::
propagate_parallel
(const tree_t &tree, const xt::xexpression<T1> &xinput)¶
-
template<typename
tree_t
, typenameT1
, typenameT2
>
autohg
::
propagate_parallel
(const tree_t &tree, const xt::xexpression<T1> &xinput, const xt::xexpression<T2> &xcondition)¶
-
template<typename
tree_t
, typenameT1
, typenameT2
>
autohg
::
propagate_sequential
(const tree_t &tree, const xt::xexpression<T1> &xinput, const xt::xexpression<T2> &xcondition)¶
-
template<typename
tree_t
, typenameT
, typenameaccumulator_t
>
autohg
::
propagate_sequential_and_accumulate
(const tree_t &tree, const xt::xexpression<T> &xinput, const accumulator_t &accumulator)¶
-
template<typename
tree_t
, typenameT1
, typenameT2
, typenameaccumulator_t
>
autohg
::
propagate_sequential_and_accumulate
(const tree_t &tree, const xt::xexpression<T1> &xinput, const accumulator_t &accumulator, const xt::xexpression<T2> &xcondition)¶
-
template<typename
graph_t
, typenametree_t
, typenameT
, typenameT1
, typenameaccumulator_t
, typenameoutput_t
= typename T::value_type>
autohg
::
accumulate_on_contours
(const graph_t &graph, const tree_t &tree, const xt::xexpression<T> &xinput, const xt::xexpression<T1> &xdepth, const accumulator_t &accumulator)¶
-
template<typename
T1
, typenameT2
>
autohg
::
project_fine_to_coarse_labelisation
(const xt::xexpression<T1> &xlabelisation_fine, const xt::xexpression<T2> &xlabelisation_coarse, size_t num_regions_fine = 0, size_t num_regions_coarse = 0)¶ Given two labelisations, a fine and a coarse one, of a same set of elements.
Find for each label (ie. region) of the fine labelisation, the label of the region in the coarse labelisation that maximises the intersection with the “fine” region.
Pre-condition: range(xlabelisation_fine) = [0..num_regions_fine[ range(xlabelisation_coarse) = [0..num_regions_coarse[
If num_regions_fine or num_regions_coarse are not provided, they will be determined as max(xlabelisation_fine) + 1 and max(xlabelisation_coarse) + 1
- Template Parameters:
T1 –
T2 –
- Parameters:
xlabelisation_fine –
num_regions_fine –
xlabelisation_coarse –
num_regions_coarse –
- Returns:
a 1d array of size num_regions_fine
-
inline auto
hg
::
project_fine_to_coarse_rag
(const region_adjacency_graph &fine_rag, const region_adjacency_graph &coarse_rag)¶ Given two region adjacency graphs, a fine and a coarse one, of a same set of elements.
Find for each region of the fine rag, the region of the coarse rag that maximises the intersection with the “fine” region.
- Parameters:
fine_rag –
coarse_rag –
- Returns:
a 1d array of size num_vertices(fine_rag.rag)
-
template<typename
graph_t
, typenameT
>
autohg
::
make_hierarchy_aligner_from_graph_cut
(const graph_t &graph, const xt::xexpression<T> &saliency_map)¶
-
template<typename
graph_t
, typenameT
>
autohg
::
make_hierarchy_aligner_from_labelisation
(const graph_t &graph, const xt::xexpression<T> &vertex_labels)¶
-
template<typename
graph_t
, typenametree_t
, typenameT
>
autohg
::
make_hierarchy_aligner_from_hierarchy
(const graph_t &graph, const tree_t &tree, const xt::xexpression<T> &altitudes)¶
-
template<typename
graph_t
, typenameT
>
autohg
::
graph_cut_2_labelisation
(const graph_t &graph, const xt::xexpression<T> &xedge_weights)¶ Labelize graph vertices according to the given graph cut.
Each edge having a non zero value in the given edge_weights are assumed to be part of the cut.
- Template Parameters:
graph_t –
T –
label_type –
- Parameters:
graph –
edge_weights –
- Returns:
-
template<typename
graph_t
, typenameT
>
autohg
::
labelisation_2_graph_cut
(const graph_t &graph, const xt::xexpression<T> &xvertex_labels)¶ Determine the graph cut that corresponds to a given labeling of the graph vertices.
The result is a weighting of the graph edges where edges with a non zero weight are part of the cut.
- Template Parameters:
graph_t –
T –
- Parameters:
graph –
xvertex_labels –
- Returns:
-
template<typename
graph_t
, typenameT
>
autohg
::
minimum_spanning_tree
(const graph_t &graph, const xt::xexpression<T> &xedge_weights)¶ Computes a minimum spanning tree of the given edge weighted graph using Kruskal’s algorithm.
The returned structure contains two elements:
the minimum spanning tree (mst)
a map (mst_edge_map) that indicates for each edge of the mst, the corresponding edge index in the input graph
If the input graph is not connected, the result is indeed a minimum spanning forest.
- Template Parameters:
graph_t – Input graph type
T – Input edge weights type
- Parameters:
graph – Input graph
xedge_weights – Input edge weights
- Returns:
a mst structure
-
template<typename
graph_t
, typenameT
>
autohg
::
subgraph_spanning
(const graph_t &graph, const xt::xexpression<T> &xedge_indices)¶ Compute a spanning subgraph of the given graph composed of the edges of the input graph indicated in the edge_indices array.
The edges of the subgraph will be in the order given in edge_indices array.
- Template Parameters:
graph_t –
T –
- Parameters:
graph – input graph
xedges_indices – list of edges of the input graph to include in the subgraph
- Returns:
a spanning subgraph
-
inline ugraph
hg
::
line_graph
(const ugraph &graph)¶ Compute the line graph of an undirected graph.
The line graph :math:
LG
of an undirected graph :math:G
is a graph such that:each vertex of :math:
LG
represents an edge of :math:G
: the :math:i
-th vertex of :math:LG
corresponds to the :math:i
-th edge of :math:G
; andtwo vertices :math:
x
and :math:y
of :math:LG
are adjacent if their corresponding edges in :math:G
share a common extremity. Formally, if :math:x
represents the edge :math:\{i, j \}
and if :math:y
represents the edge :math:\{k, j \}
, then the edge :math:\{x, y\}
belongs to :math:LG
if :math:\{i, j \} \\cap \{k, j \} \\neq \emptyset
.
The line graph is also known as: the covering graph, the derivative, the edge-to-vertex dual, the conjugate, the representative graph, the edge graph, the interchange graph, the adjoint graph, or the derived graph.
- Parameters:
graph –
- Returns:
-
template<typename
graph_t
>
ugraphhg
::
line_graph
(const graph_t &graph)¶ See description of the function above.
- Template Parameters:
graph_t –
- Parameters:
graph –
- Returns:
-
template<typename
result_value_t
= double, typenamegraph_t
>
autohg
::
weight_graph
(const graph_t &graph, const std::function<result_value_t(typename graph_t::vertex_descriptor, typename graph_t::vertex_descriptor)> &fun)¶ Compute edge-weights of a graph based on a weighting function.
A weighting function is a function that associates a weights to a pair of vertices.
- Template Parameters:
graph_t –
result_value_t –
- Parameters:
graph –
fun –
- Returns:
an array of weights
-
template<typename
result_value_t
= double, typenamepromoted_type
= double, typenamegraph_t
, typenameT
>
autohg
::
weight_graph
(const graph_t &graph, const xt::xexpression<T> &xvertex_weights, weight_functions weight)¶ Compute edge-weights of a graph based from the vertex-weights and a predefined weighting function (see weight_functions enum).
Each edge is weighted with a combination of its extremities weights.
- Template Parameters:
result_value_t – The value type of the result
promoted_type – The value type used for internal computation
graph_t –
T –
- Parameters:
graph –
xvertex_weights –
weight –
- Returns:
-
template<typename
value_t
>
decltype(auto)hg
::
make_horizontal_cut_nodes
(array_1d<index_t> &&nodes, value_t altitude)¶
-
template<typename
tree_t
, typenameT
>
decltype(auto)hg
::
make_horizontal_cut_explorer
(tree_t &&tree, T &&altitudes)¶
-
template<typename
graph_t
, typenameT
>
autohg
::
make_region_adjacency_graph_from_labelisation
(const graph_t &graph, const xt::xexpression<T> &xvertex_labels)¶ Construct a region adjacency graph from a vertex labeled graph in linear time.
- Template Parameters:
graph_t –
T –
- Parameters:
graph –
xvertex_labels –
- Returns:
see struct region_adjacency_graph
-
template<typename
graph_t
, typenameT
>
autohg
::
make_region_adjacency_graph_from_graph_cut
(const graph_t &graph, const xt::xexpression<T> &xedge_weights)¶ Construct a region adjacency graph from a graph cut in linear time.
Any edge with weight different from 0 belongs to the cut.
TODO factorize method with make_region_adjacency_graph_from_labelisation
- Template Parameters:
graph_t –
T –
- Parameters:
graph –
xedge_weights –
- Returns:
see struct region_adjacency_graph
-
template<typename
T
>
autohg
::
rag_back_project_weights
(const array_1d<index_t> &rag_map, const xt::xexpression<T> &xrag_weights)¶ Projects weights on the rag (vertices or edges) onto the original graph.
- Template Parameters:
T –
- Parameters:
rag_map – rag vertex_map or rag edge_map (see struct region_adjacency_graph)
xrag_weights – node or edge weights of the rag (depending of the provided rag_map)
- Returns:
original graph (vertices or edges) weights
-
template<typename
T
, typenameaccumulator_t
, typenameoutput_t
= typename T::value_type>
autohg
::
rag_accumulate
(const array_1d<index_t> &rag_map, const xt::xexpression<T> &xweights, const accumulator_t &accumulator)¶ Accumulate original graph (vertices or edges) weights onto the rag (vertices or edges)
- Template Parameters:
T –
accumulator_t –
output_t –
- Parameters:
rag_map – rag vertex_map or rag edge_map (see struct region_adjacency_graph)
xweights – node or edge weights of the original graph (depending of the provided rag_map)
accumulator –
- Returns:
rag (vertices or edges) weights
-
template<typename
tree_t
, typenameT1
, typenameT2
>
autohg
::
reconstruct_leaf_data
(const tree_t &tree, const xt::xexpression<T1> &altitudes, const xt::xexpression<T2> &deleted_nodes)¶ Each leaf of the tree takes the weight of its closest non deleted ancestor.
- Template Parameters:
tree_t –
T1 –
T2 –
- Parameters:
tree –
altitudes –
deleted_nodes –
- Returns:
-
template<typename
tree_t
, typenameT
, typenamevalue_t
>
autohg
::
labelisation_horizontal_cut_from_threshold
(const tree_t &tree, const xt::xexpression<T> &xaltitudes, const value_t threshold)¶ Labelize tree leaves according to an horizontal cut in the tree.
Two leaves are in the same region (ie. have the same label) if the altitude of their lowest common ancestor is smaller or equal than the specified threshold.
The label of a leave l is equal to the index of the smallest node containing l and whose altitude is strictly greater than the specified threshold.
- Template Parameters:
tree_t –
T –
value_t –
- Parameters:
tree –
xaltitudes –
threshold –
- Returns:
-
template<typename
tree_t
, typenameT
>
autohg
::
labelisation_hierarchy_supervertices
(const tree_t &tree, const xt::xexpression<T> &xaltitudes)¶ Labelize the tree leaves into supervertices.
Two leaves are in the same supervertex if they have a common ancestor of altitude 0.
This functions guaranties that the labels are in the range [0, num_supervertices-1].
- Template Parameters:
tree_t –
T –
- Parameters:
tree –
xaltitudes –
- Returns:
-
template<typename
tree_t
, typenameT
>
autohg
::
supervertices_hierarchy
(const tree_t &tree, const xt::xexpression<T> &xaltitudes)¶ Extract the supervertices associated to the given tree and create the equivalent tree on this supervertices.
Two leaves are in the same supervertex if they have a common ancestor of altitude 0.
The equivalent tree is obtained by removing all nodes of the given tree which does not contain any of the supervertices. Its leaves are thus the supervertices.
Also returns an array that maps any node index i of the new tree, to the index of this node in the original tree.
- Template Parameters:
tree_t –
T –
- Parameters:
tree –
xaltitudes –
- Returns:
-
template<typename
tree1_t
, typenametree2_t
>
boolhg
::
test_tree_isomorphism
(const tree1_t &t1, const tree2_t &t2)¶ Test if 2 trees are isomorph assuming that they share the same leaves.
By this definition t1 is isomorph to t2 if there exist a bijection f from vertices(t1) to vertices(t2) such that: 1) for any leaf node n of t1, f(n) = n and n 2) for any node n of t1, f(t1.parent(n)) = t2.parent(f(n))
Note that the root r node of a tree t is defined by t.parent(r) = r, thus 2) becomes for the root node r1 of t1, f(r1) = t2.parent(f(r1)), i.e. f(r1) is the root node of t2
- Template Parameters:
tree1_t –
tree2_t –
- Parameters:
t1 –
t2 –
- Returns:
-
template<typename
tree_t
, typenameT1
, typenameT2
>
autohg
::
binary_labelisation_from_markers
(const tree_t &tree, const xt::xexpression<T1> &xobject_marker, const xt::xexpression<T2> &xbackground_marker)¶ Given two binary markers o (object) and b (background) (given by their indicator functions) on the leaves of a tree t, the corresponding binary labelization of the leaves of t is defined as the union of all the nodes intersecting o but not b.
final_object = union {R in T | R cap o neq emptyset and R cap b = emptyset}
- Template Parameters:
tree_t – tree type
T1 – xtensor type, value_type must be castable to bool
T2 – xtensor type, value_type must be castable to bool
- Parameters:
tree – input tree
xobject_marker – indicator function of the object marker
xbackground_marker – indicator function of the background marker
- Returns:
indicator function of the final_object
-
template<typename
tree_t
, typenameT
>
autohg
::
sort_hierarchy_with_altitudes
(const tree_t &tree, const xt::xexpression<T> &xaltitudes)¶ Sort the nodes of a tree according to their altitudes.
The altitudes must be increasing, i.e. for any nodes i, j such that j is an ancestor of j, then altitudes[i] <= altitudes[j].
The result is a new tree and a node map, isomorph to the input tree such that for any nodes i and j, i<j => altitudes[node_map[i]] <= altitudes[node_map[j]]
The latter condition is stronger than the original condition on the altitudes as j is an ancestor of i implies i<j while the converse is not true.
Note that the altitudes of the new tree can be obtained with:
auto res = sort_hierarchy_with_altitudes(tree, altitudes); auto & new_tree = res.tree; auto new_altitudes = xt::index_view(altitudes, res.node_map);
- Template Parameters:
tree_t –
T –
- Parameters:
tree –
xaltitudes –
- Returns:
-
template<typename
tree_t
>
autohg
::
sub_tree
(const tree_t &tree, index_t root)¶ Extract the sub tree rooted in the given node from the given tree.
The result is a new tree :math:
st
and a node map :math:nm
such thatthe node map associates each node of the sub tree :math:
st
to its corresponding node in the original treethe order of the nodes of the original tree is preserved in the sub tree: for any vertices :math:
x
and :math:y
of :math:st
such that :math:x < y
then :math:nm[x] < nm[y]
:Complexity:
The tree is constructed in linearithmic time :math:
\mathcal{O}(n\log(n))
with :math:n
the number of vertices in the sub tree.- Template Parameters:
tree_t –
- Parameters:
tree –
root –
- Returns:
-
template<typename
tree_type
, typenameT
, typenameaccumulator_type
= hg::accumulator_sum>
autohg
::
labelisation_optimal_cut_from_energy
(const tree_type &tree, const xt::xexpression<T> &xenergy_attribute, const accumulator_type accumulator = hg::accumulator_sum())¶ Computes the labelisation of the input tree leaves corresponding to the optimal cut according to the given energy attribute.
Given a node i, the value energy_attribute(i) represents the energy fo the partial partition composed of the single region i. Given a node i, the energy of the partial partition composed of the children of i is given by accumulator(energy_attribute(children(i))).
This function computes the partition (ie. a set of node forming a cut of the tree) that has a minimal energy according to the definition above.
The algorithm used is based on dynamic programming and runs in linear time w.r.t. to the number of nodes in the tree.
See:
Laurent Guigues, Jean Pierre Cocquerez, H ervé Le Men. Scale-sets Image Analysis. International Journal of Computer Vision, Springer Verlag, 2006, 68 (3), pp.289-317
and
Bangalore Ravi Kiran, Jean Serra. Global-local optimizations by hierarchical cuts and climbing energies. Pattern Recognition Letters, Elsevier, 2014, 47 (1), pp.12-24.
- Template Parameters:
tree_type – input tree type
T – energy attribute type
accumulator_type – accumulator type
- Parameters:
tree – input tree
xenergy_attribute – 1d array of energy attribute for the input tree
accumulator – accumulator used to define how children energies are combined in order to obtain the energy of the corresponding partial partition
- Returns:
a 1d integer array with num_leaves(tree) elements representing the minimal energy partition
-
template<typename
tree_type
, typenameT
>
autohg
::
hierarchy_to_optimal_energy_cut_hierarchy
(const tree_type &tree, const xt::xexpression<T> &xdata_fidelity_attribute, const xt::xexpression<T> &xregularization_attribute, const int approximation_piecewise_linear_function = 10)¶ Transforms the given hierarchy into its optimal energy cut hierarchy for the given energy terms.
In the optimal energy cut hierarchy, any horizontal cut corresponds to an optimal energy cut in the original hierarchy.
Each node i of the tree is associated to a data fidelity energy D(i) and a regularization energy R(i). The algorithm construct a new hierarchy with associated altitudes such that the horizontal cut of level lambda is the optimal cut for the energy attribute D + lambda * R of the input tree (see function labelisation_optimal_cut_from_energy). In other words, the horizontal cut of level lambda in the result is the cut of the input composed of the nodes N such that sum_{r in N} D(r) + lambda * R(r) is minimal.
PRECONDITION: the regularization energy R must be sub additive: for each node i: R(i) <= sum_{c in children(i)} R(c)
The algorithm runs in linear time O(n)
See:
Laurent Guigues, Jean Pierre Cocquerez, Hervé Le Men. Scale-sets Image Analysis. International Journal of Computer Vision, Springer Verlag, 2006, 68 (3), pp.289-317
- Template Parameters:
tree_type –
T –
- Parameters:
tree – Input tree
xdata_fidelity_attribute – Data fidelity energy (1d array)
xregularization_attribute – Regularization energy (1d array)
approximation_piecewise_linear_function – Maximum number of pieces used in the approximated piecewise linear model for the energy.
- Returns:
-
template<typename
graph_t
, typenameT1
, typenameT2
, typenameT3
, typenameT4
, typenameT5
>
autohg
::
binary_partition_tree_MumfordShah_energy
(const graph_t &graph, const xt::xexpression<T1> &xvertex_perimeter, const xt::xexpression<T2> &xvertex_area, const xt::xexpression<T3> &xvertex_values, const xt::xexpression<T4> &xsquared_vertex_values, const xt::xexpression<T5> &xedge_length)¶ Compute the binary partition tree, i.e.
the agglomerative clustering, according to the Mumford-Shah energy with a constant piecewise model.
The distance between two regions is equal to the apparition scale of the merged region.
See:
Laurent Guigues, Jean Pierre Cocquerez, Hervé Le Men. Scale-sets Image Analysis. International Journal of Computer Vision, Springer Verlag, 2006, 68 (3), pp.289-317
- Template Parameters:
graph_t –
T1 –
T2 –
T3 –
T4 –
T5 –
- Parameters:
graph – Input graph
xvertex_perimeter – Perimeter of each vertex of the input graph
xvertex_area – Area of each vertex of the input graph
xvertex_values – Sum of values inside the region represented by each vertex of the input graph.
xsquared_vertex_values – Sum of the squared values inside the region represented by each vertex of the input graph.
xedge_length – Length of the frontier represented by each edge.
- Returns:
-
template<typename
tree_iterator
>
autohg
::
tree_fusion_depth_map
(const tree_iterator first, const tree_iterator last)¶ The result is a directed acyclic graph with a single root (corresponding to the roots of the input trees).
The depth of a node in this graph is defined as the length of the longest path from the root this node.
This function returns the depth of the leaves of this graph (which are the same as the leaves of the input trees).
- Template Parameters:
tree_iterator – Iterator on tree pointers
- Parameters:
first –
last –
- Returns:
an 1d array of size num_leaves(**first)
-
template<typename
range_tree_t
>
autohg
::
tree_fusion_depth_map
(const range_tree_t &range)¶
-
template<typename
tree_t
, typenameT
, typenameTw
>
autohg
::
tree_monotonic_regression
(const tree_t &tree, const xt::xexpression<T> &xaltitudes, const xt::xexpression<Tw> &xweights, const std::string &mode)¶ Monotonic regression on the given tree altitudes.
Computes new altitudes
naltitudes
that are close to the given :attr:altitudes
and that are increasing for the given :attr:tree
: i.e. for any nodes :math:i, j
such that :math:j
is an ancestor of :math:i
, then :math:naltitudes[i] \leq naltitudes[j]
.The definition of close depends of the value of :attr:
mode
:If :attr:
mode
is equal to"min"
thennaltitudes
is the largest increasing function below :attr:altitudes
.If :attr:
mode
is equal to"max"
thennaltitudes
is the smallest increasing function above :attr:altitudes
.If :attr:
mode
is equal to"least_square"
thennaltitudes
minizes the following minization problem:
.. math::
naltitudes = \\arg \min_x \sum_i (weights[i] * (altitudes[i] - x[i])^2)
such that
naltitudes
is increasing for :attr:tree
.:Complexity:
With :math:
n
the number of nodes in the :attr:tree
:For the modes
"min"
and"max"
, the runtime complexity is linear :math:\mathcal{O}(n)
.For the mode
"least_square"
, the runtime complexity is linearithmic :math:\mathcal{O}(n\log(n))
and the space complexity is linear :math:\mathcal{O}(n)
. The algorithm used is described in:P. Pardalos and G. Xue ‘Algorithms for a Class of Isotonic Regression Problems. https://link.springer.com/article/10.1007/PL00009258`_ Algorithmica (1999) 23: 211. doi:10.1007/PL00009258
- Template Parameters:
tree_t –
T –
Tw –
- Parameters:
tree – input tree
xaltitudes – tree node altitudes
xweights – tree node weights
mode – “min”, “max”, or ” least_square”
- Returns:
-
template<typename
tree_t
, typenameT
>
autohg
::
tree_monotonic_regression
(const tree_t &tree, const xt::xexpression<T> &xaltitudes, const std::string &mode)¶
-
template<typename
graph_t
, typenameT
>
autohg
::
labelisation_watershed
(const graph_t &graph, const xt::xexpression<T> &xedge_weights)¶ Linear time watershed cut algorithm.
Jean Cousty, Gilles Bertrand, Laurent Najman, Michel Couprie. Watershed Cuts: Minimum Spanning Forests and the Drop of Water Principle. IEEE Transactions on Pattern Analysis and Machine Intelligence, Institute of Electrical and Electronics Engineers, 2009, 31 (8), pp.1362-1374.
- Template Parameters:
graph_t –
T –
- Parameters:
graph –
xedge_weights –
- Returns:
array of labels on graph vertices, numbered from 1 to n with n the number of minima
-
template<typename
graph_t
, typenameT1
, typenameT2
>
autohg
::
labelisation_seeded_watershed
(const graph_t &graph, const xt::xexpression<T1> &xedge_weights, const xt::xexpression<T2> &xvertex_seeds, const typename T2::value_type background_label = 0)¶
-
template<typename
tree_t
, typenameT
>
autohg
::
dendrogram_purity
(const tree_t &tree, const xt::xexpression<T> &xleaf_labels)¶ Weighted average of the purity of each node of the tree with respect to a ground truth labelization of the tree leaves.
Let :math:
T
be a tree with leaves :math:V=\{1, \ldots, n\}
. Let :math:C=\{C_1, \ldots, C_K\}
be a partition of :math:V
into :math:k
(label) sets.The purity of a subset :math:
X
of :math:V
with respect to class :math:C_\ell\in C
is the fraction of elements of :math:X
that belongs to class :math:C_\ell
:.. math::
pur(X, C_\ell) = \\frac{| X \cap C_\ell |}{| X |}.
The purity of :math:
T
is the defined as:.. math::
pur(T) = \\frac{1}{Z}\sum_{k=1}^{K}\sum_{x,y\in C_k, x\\neq y} pur(lca_T(x,y), C_k)
with :math:
Z=| \{\{x,y\} \subseteq V \mid x\\neq y, \exists k, \{x,y\}\subseteq C_k\} |
.:See:
Heller, Katherine A., and Zoubin Ghahramani. "`Bayesian hierarchical clustering <https://www2.stat.duke.edu/~kheller/bhcnew.pdf>`_ ." Proc. ICML. ACM, 2005.
:Complexity:
The dendrogram purity is computed in :math:
\mathcal{O}(N\\times K \\times C^2)
with :math:N
the number of nodes in the tree, :math:K
the number of classes, and :math:C
the maximal number of children of a node in the tree.Let :maht:
T
be a tree with leaves :maht:V=\{1, \ldots, n\}
. Let :maht:C=\{C_1, \ldots, C_K\}
be a partition of :maht:V
into :maht:k
(label) sets.The purity of a subset :maht:
X
of :maht:V
with respect to class :maht:C_\ell\in C
is the fraction of elements of :maht:X
that belongs to class :maht:C_\ell
:.. math::
pur(X, C_\ell) = \\frac{| X \cap C_\ell |}{| X |}.
The purity of :maht:
T
is the defined as:.. math::
pur(T) = \\frac{1}{Z}\sum_{k=1}^{K}\sum_{x,y\in C_k, x\neq y} pur(lca_T(x,y), C_k)
with :math:
Z=| \{\{x,y\} \subseteq V \mid x\neq y, \exists k, \{x,y\}\subseteq C_k\} |
.:See:
Heller, Katherine A., and Zoubin Ghahramani. "Bayesian hierarchical clustering." Proc. ICML. ACM, 2005.
:Complexity:
The dendrogram purity is computed in :math:
\mathcal{O}(N\times K \times C^2)
with :math:N
the number of nodes in the tree, :math:K
the number of classes, and :math:C
the maximal number of children of a node in the tree.- Template Parameters:
tree_t –
T –
tree_t –
T –
- Parameters:
tree – input tree
xleaf_labels – must be a 1d array with values in [0, max_label]
tree – input tree
xleaf_labels – must be a 1d array with values in [0, max_label]
- Returns:
a score between 0 and 1 (higher is better)
- Returns:
a score between 0 and 1 (higher is better)
-
template<typename
tree_t
, typenameT1
, typenameT2
, typenamescorer_t
>
autohg
::
assess_fragmentation_horizontal_cut
(const tree_t &tree, const xt::xexpression<T1> &xaltitudes, const xt::xexpression<T2> &xground_truth, const scorer_t &partition_scorer, const array_1d<index_t> &vertex_map = {}, size_t max_regions = 200)¶
-
template<typename
value_type
= index_t, typenameT1
, typenameT2
>
autohg
::
card_intersections
(const xt::xexpression<T1> &xcandidate, const xt::xexpression<T2> &xground_truths)¶
-
template<typename
T
, typenamescorer_t
>
autohg
::
assess_partition
(const std::vector<T> &card_intersections, const scorer_t &scorer)¶
-
template<typename
T1
, typenameT2
, typenamescorer_t
>
autohg
::
assess_partition
(const xt::xexpression<T1> &xcandidate, const xt::xexpression<T2> &xground_truths, const scorer_t &scorer)¶
-
template<typename
tree_t
, typenameT
>
autohg
::
attribute_area
(const tree_t &tree, const xt::xexpression<T> &xleaf_area)¶ The area of a node n of the tree t is equal to the sum of the area of the leaves in the subtree rooted in n.
area(n) = sum_{l in leaves(t), l is a descendant of n} area(n)
- Template Parameters:
tree_t – tree type
T – xexpression derived type of xleaf_area
- Parameters:
tree – input tree
xleaf_area – area of the leaves of the input tree
- Returns:
an array with the area of each node of the tree
-
template<typename
tree_t
>
autohg
::
attribute_area
(const tree_t &tree)¶ The area of a node n of the tree t is equal to the number of leaves in the subtree rooted in n.
area(n) = |{l in leaves(t), l is a descendant of n}|
- Template Parameters:
tree_t – tree type
- Parameters:
tree – input tree
- Returns:
an array with the area of each node of the tree
-
template<typename
tree_t
, typenameT1
, typenameT2
>
autohg
::
attribute_volume
(const tree_t &tree, const xt::xexpression<T1> &xnode_altitude, const xt::xexpression<T2> &xnode_area)¶ The volume of a node n of the tree t is defined recursively as: volume(n) = abs(altitude(n) - altitude(parent(n)) * area(n) + sum_{c in children(n, t)} volume(c)
- Template Parameters:
tree_t – tree type
T1 – xexpression derived type of xnode_altitude
T2 – xexpression derived type of xnode_area
- Parameters:
tree – input tree
xnode_altitude – altitude of the nodes of the input tree
xnode_area – area of the nodes of the input tree
- Returns:
an array with the volume of each node of the tree
-
template<typename
tree_t
>
autohg
::
attribute_depth
(const tree_t &tree)¶ The depth of a node n of the tree t is equal to the number of ancestors of n.
- Template Parameters:
tree_t – tree type
- Parameters:
tree – input tree
- Returns:
an array with the depth of each node of the tree
-
template<typename
tree_t
, typenameT
>
autohg
::
attribute_height
(const tree_t &tree, const xt::xexpression<T> &xaltitudes, bool increasing_altitudes)¶ In a tree :math:
t
, given that the altitudes of the nodes vary monotically from the leaves to the root, the height of a node :math:n
of :math:t
is equal to the difference between the altitude of the parent of :math:n
and the altitude of the deepest non-leaf node in the subtree of :math:t
rooted in :math:n
.If :attr:
increasing_altitude
is true, this means that altitudes are increasing from the leaves to the root (ie. for any node :math:n
, :math:altitudes(n) \leq altitudes(parent(n))
. Else, if :attr:increasing_altitude
is false, this means that altitudes are decreasing from the leaves to the root (ie. for any node :math:n
, :math:altitude(n) \geq altitude(parent(n))
.- Template Parameters:
tree_t – tree type
T – xexpression derived type of xnode_altitude
- Parameters:
tree – input tree
xnode_altitude – altitude of the nodes of the input tree
increasing_altitudes – must be true if altitude is increasing, false if it is decreasing
- Returns:
an array with the height of each node of the tree
-
template<typename
tree_t
, typenameT
>
autohg
::
attribute_extrema
(const tree_t &tree, const xt::xexpression<T> &xaltitudes)¶ Identify nodes in a hierarchy that represent extrema.
An extremum of the hierarchy :math:
T
with altitudes :math:alt
is a node :math:n
of :math:T
such that the altitude of any non leaf node included in :math:n
is equal to the altitude of :math:n
and the altitude of the parent of :math:n
is different from the altitude of :math:n
.The result is a boolean array such that :math:
result(n)=true
if the node :math:n
is an extremum.- Template Parameters:
tree_t –
T –
- Parameters:
tree – Input tree
xaltitudes – Tree node altitudes
- Returns:
-
template<typename
tree_t
, typenameT1
, typenameT2
>
autohg
::
attribute_extinction_value
(const tree_t &tree, const xt::xexpression<T1> &xaltitudes, const xt::xexpression<T2> &xattribute, bool increasing_altitudes)¶ The extinction value of a node :math:
n
of the input tree :math:t
with increasing altitudes :math:alt
for the increasing attribute :math:att
is the equal to the threshold :math:k
such that the node :math:n
is still in an minima of :math:t
when all nodes having an attribute value smaller than :math:k
are removed.Formally, let :math:
\{M_i\}
be the set of minima of :math:(t, altitudes)
. Let :math:prec
be a total ordering of :math:\{M_i\}
such that :math:M_i \prec M_j \Rightarrow alt(M_i) \leq alt(M_j)
. Let :math:r(M_i)
be the smallest node of :math:t
containing :math:M_i
and another minima :mathM_j
such that :math:M_j \prec M_i
. The extinction value of :mathM_i
is then defined as :math:alt(r(M_i)) - alt(M_i)
.Extinction values of minima are then extended to other nodes in the tree with the following rules:
the extinction value of a non-leaf node :math:
n
which is not a minimum is defined as the largest extinction values among all the minima contained in :math:n
(and 0 if :math:n
does not contain any minima); andthe extinction value of a leaf node :math:
n
belonging to a minima :math:M_i
is equal to the extinction value of :math:M_i
. I :math:n
does not belong to any minima its extinction value is 0.
If :attr:
increasing_altitude
is true, this means that altitudes are increasing from the leaves to the root (ie. for any node :math:n
, :math:altitudes(n) \leq altitudes(parent(n))
. Else, if :attr:increasing_altitude
is false, this means that altitudes are decreasing from the leaves to the root (ie. for any node :math:n
, :math:altitude(n) \geq altitude(parent(n))
: you should then replace minima by maxima in the description above.- Template Parameters:
tree_t – tree type
T – xexpression derived type of xaltitude
- Parameters:
tree – input tree
xaltitudes – altitude of the nodes of the input tree
xattribute – attribute used for filtering
increasing_altitudes – must be true if altitude is increasing, false if it is decreasing
- Returns:
an array with the dynamics of each node of the tree
-
template<typename
tree_t
, typenameT
>
autohg
::
attribute_dynamics
(const tree_t &tree, const xt::xexpression<T> &xaltitudes, bool increasing_altitudes)¶ Given a node :math:
n
of the tree :math:t
, the dynamics of :math:n
is the difference between the altitude of the deepest minima of the subtree rooted in :math:n
and the altitude of the closest ancestor of :math:n
that has a deeper minima in its subtree.If no such ancestor exists then, the dynamics of :math:
n
is equal to the difference between the altitude of the highest node of the tree (the root) and the depth of the deepest minima.The dynamics is the extinction values for the attribute height.
If :attr:
increasing_altitude
is true, this means that altitudes are increasing from the leaves to the root (ie. for any node :math:n
, :math:altitudes(n) \leq altitudes(parent(n))
. Else, if :attr:increasing_altitude
is false, this means that altitudes are decreasing from the leaves to the root (ie. for any node :math:n
, :math:altitude(n) \geq altitude(parent(n))
.- Template Parameters:
tree_t – tree type
T – xexpression derived type of xaltitudes
- Parameters:
tree – input tree
xaltitudes – altitude of the nodes of the input tree
increasing_altitudes – must be true if altitude is increasing, false if it is decreasing
- Returns:
an array with the dynamics of each node of the tree
-
template<typename
tree_t
>
autohg
::
attribute_sibling
(const tree_t &tree, index_t skip = 1)¶ For each node
n
which is thek
-th child of its parent nodep
amongN
children, the attribute sibling ofn
is the index of the(k + skip) % N
-th child ofp
.The sibling of the root node is itself.
The sibling attribute enables to easily emulates a (doubly) linked list among brothers.
In a binary tree, the sibling attribute of a node is effectively its only brother (with
skip
equals to 1).- Template Parameters:
tree_t –
- Parameters:
tree – Input tree
skip – Number of skipped element in the children list (including yourself)
- Returns:
an array with the sibling index of each node of the tree
-
template<typename
tree_t
, typenamegraph_t
, typenameT1
, typenameT2
>
autohg
::
attribute_contour_length_component_tree
(const tree_t &tree, const graph_t &base_graph, const xt::xexpression<T1> &xvertex_perimeter, const xt::xexpression<T2> &xedge_length)¶ Computes the contour length (perimeter) of each node of the input component tree.
Warning: does not work for tree of shapes left in original space (the problem is that two children of a node may become adjacent when the interpolated pixels are removed).
- Template Parameters:
tree_t –
graph_t –
T1 –
T2 –
- Parameters:
tree – input tree
base_graph – graph on the leaves of tree
xvertex_perimeter – perimeter of each vertex of the base graph
xedge_length – length of each edge of the base graph (length of the frontier between the two adjacent vertices)
- Returns:
-
template<typename
tree_t
>
autohg
::
attribute_child_number
(const tree_t &tree)¶ Given a node :math:
n
whose parent is :math:p
, the attribute value of :math:n
is the rank of :math:n
in the list of children of :math:p
.In other :math:
attribute(n)=i
means that :math:n
is the :math:i
-th child of :math:p
.The root of the tree, who has no parent, take the value -1.
- Template Parameters:
tree_t –
- Parameters:
tree –
- Returns:
-
template<typename
tree_t
>
autohg
::
attribute_smallest_enclosing_shape
(const tree_t &t1, const tree_t &t2)¶ Given two trees :math:
t_1
and :math:t_2
defined over the same domain, ie sharing the same set of leaves.For each node :math:
n
of :math:t1
, computes the index of the smallest node of :math:t2
containing :math:n
.- Template Parameters:
tree_t –
- Parameters:
t1 –
t2 –
- Returns:
-
template<typename
tree_t
, typenameT
, typenamevalue_type
= typename T::value_type>
autohg
::
attribute_children_pair_sum_product
(const tree_t &tree, const xt::xexpression<T> &xnode_weights)¶ Given a tree :math:
T
with node weights :math:w
: the children pair sum product for a node :math:n
sums for every pairs :math:(c_i, c_j)
of children of :math:n
, the product of the node weights of :math:c_i
and :math:c_j
.Formally:
.. math::
res(n) = \sum_{i=0}^{i<numc(n)} \sum_{j=0}^{j<i} w(child(i, n)) * w(child(j, n))
where :math:
numc(n)
is the number of children of :math:n
and :math:child(i, n)
is the :math:i
-th child of the node :math:n
.The result is thus an array of the same shape of :attr:
node_weights
- Template Parameters:
tree_t –
T –
value_type –
- Parameters:
tree –
xnode_weights –
- Returns:
-
template<typename
graph_t
>
autohg
::
source
(const std::pair<typename graph::graph_traits<graph_t>::vertex_descriptor, typename graph::graph_traits<graph_t>::vertex_descriptor> &e, const graph_t&)¶ Source vertex of an edge.
- Template Parameters:
graph_t –
- Parameters:
e –
- Returns:
-
template<typename
graph_t
>
autohg
::
sources
(const graph_t &t)¶ Source vertex of all edges oin the given graph.
- Template Parameters:
graph_t –
- Parameters:
t –
- Returns:
1d expression with num_edges(t) element
-
template<typename
graph_t
>
autohg
::
target
(const std::pair<typename graph::graph_traits<graph_t>::vertex_descriptor, typename graph::graph_traits<graph_t>::vertex_descriptor> &e, const graph_t&)¶ Target vertex of an edge.
- Template Parameters:
graph_t –
- Parameters:
e –
- Returns:
-
template<typename
graph_t
>
autohg
::
targets
(const graph_t &t)¶ Target vertex of all edges oin the given graph.
- Template Parameters:
graph_t –
- Parameters:
t –
- Returns:
1d expression with num_edges(t) element
-
template<typename
graph_t
>
autohg
::
vertex_iterator
(const graph_t &g)¶ Range over all vertices of the given graph.
- Template Parameters:
graph_t –
- Parameters:
g –
- Returns:
-
template<typename
graph_t
>
autohg
::
edge_iterator
(const graph_t &g)¶ Range over all edges of the given graph.
- Template Parameters:
graph_t –
- Parameters:
g –
- Returns:
-
template<typename
graph_t
>
autohg
::
out_edge_iterator
(typename graph::graph_traits<graph_t>::vertex_descriptor v, const graph_t &g)¶ Range over all edges whose source is the given vertex in the given graph.
- Template Parameters:
graph_t –
- Parameters:
v –
g –
- Returns:
-
template<typename
graph_t
>
autohg
::
in_edge_iterator
(typename graph::graph_traits<graph_t>::vertex_descriptor v, const graph_t &g)¶ Range over all edges whose target is the given vertex in the given graph.
- Template Parameters:
graph_t –
- Parameters:
v –
g –
- Returns:
-
template<typename
graph_t
>
autohg
::
adjacent_vertex_iterator
(typename graph::graph_traits<graph_t>::vertex_descriptor v, const graph_t &g)¶ Range over all vertices adjacent to the given vertex.
- Template Parameters:
graph_t –
- Parameters:
v –
g –
- Returns:
-
template<typename
graph_t
>
autohg
::
children_iterator
(typename graph_t::vertex_descriptor v, const graph_t &g)¶ Range over the children vertices of the given node in the given tree.
- Template Parameters:
graph_t –
- Parameters:
v –
g –
- Returns:
-
template<typename
graph_t
>
autohg
::
ancestors_iterator
(typename graph_t::vertex_descriptor v, const graph_t &g)¶ Range over the ancestors of v in topological order (starting from v included)
- Template Parameters:
graph_t –
- Parameters:
v –
g –
- Returns:
-
template<typename
T
, typenamegraph_t
>
autohg
::
degree
(const xt::xexpression<T> &xindex, const graph_t &g)¶ Degrees of all the given vertices in the given graph.
- Template Parameters:
T – type of indices (must be integral, preferably index_t)
graph_t –
- Parameters:
xindex – array of vertex indices
g –
- Returns:
array of the same size as xindex containing the degree of each vertex indicated in xindex
-
template<typename
T
, typenamegraph_t
>
autohg
::
in_degree
(const xt::xexpression<T> &xindex, const graph_t &g)¶ In-degrees of all the given vertices in the given graph.
- Template Parameters:
T – type of indices (must be integral, preferably index_t)
graph_t –
- Parameters:
xindex – array of vertex indices
g –
- Returns:
array of the same size as xindex containing the in-degree of each vertex indicated in xindex
-
template<typename
T
, typenamegraph_t
>
autohg
::
out_degree
(const xt::xexpression<T> &xindex, const graph_t &g)¶ Out-degrees of all the given vertices in the given graph.
- Template Parameters:
T – type of indices (must be integral, preferably index_t)
graph_t –
- Parameters:
xindex – array of vertex indices
g –
- Returns:
array of the same size as xindex containing the out-degree of each vertex indicated in xindex
-
template<typename
T
, typenamegraph_t
>
voidhg
::
add_edges
(const xt::xexpression<T> &xsources, const xt::xexpression<T> &xtargets, graph_t &g)¶ Add all edges given as a pair of arrays (sources, targets) to the graph.
- Template Parameters:
T – xexpression type
graph_t – Mutable graph type
- Parameters:
xsources – Must be a 1d array of integral values
xtargets – Must have the same shape as xsources
g – A mutable graph
-
template<typename
output_graph_type
= ugraph, typenameT
>
output_graph_typehg
::
copy_graph
(const T &graph)¶ Create a new graph as a copy of the given graph.
- Template Parameters:
T – input graph type
output_graph_type – return type (default = ugraph)
- Parameters:
graph –
- Returns:
-
template<typename
output_graph_type
>
output_graph_typehg
::
copy_graph
(const ugraph &graph)¶
-
template<typename
graph_t
>
autohg
::
other_vertex
(const typename graph::graph_traits<graph_t>::edge_descriptor &edge, typename graph::graph_traits<graph_t>::vertex_descriptor vertex, const graph_t &graph)¶ Given an edge and one of the two extremities of this edge, return the other extremity (if the source is given it returns the target and vice versa).
- Template Parameters:
graph_t –
- Parameters:
edge –
vertex –
graph –
- Returns:
-
template<typename
undirected_graph
, typenameT
, typenamevalue_type
= typename T::value_type>
autohg
::
undirected_graph_2_adjacency_matrix
(const undirected_graph &graph, const xt::xexpression<T> &xedge_weights, const value_type &non_edge_value = 0)¶ Create an adjacency matrix from an undirected edge-weighted graph (the result is thus symmetric).
As the given graph is not necessarily complete, non-existing edges will receive the value
non_edge_value
in the adjacency matrix.- Template Parameters:
undirected_graph –
T –
value_type –
- Parameters:
graph – Input undirected graph
xedge_weights – Input edge-weights
non_edge_value – Value used to represent non existing edges
- Returns:
A 2d square array
-
template<typename
T
, typenamevalue_type
= typename T::value_type>
autohg
::
adjacency_matrix_2_undirected_graph
(const xt::xexpression<T> &xadjacency_matrix, const value_type &non_edge_value = 0)¶ Creates an undirected edge-weighted graph from an adjacency matrix.
Adjacency matrix entries which are equal to
non_edge_value
are not considered to be part of the graph.- Template Parameters:
T –
value_type –
- Parameters:
xadjacency_matrix – Input adjacency matrix
non_edge_value – Value used to represent non existing edges
- Returns:
a pair of types (ugraph, array_1d) representing the graph and its edge-weights
-
template<typename
graph_t
, typenameweighter
, typenameT
>
autohg
::
binary_partition_tree
(const graph_t &graph, const xt::xexpression<T> &xedge_weights, weighter weight_function)¶ Compute the binary partition tree of the graph.
At each step: 1 - the algorithm finds the edge of smallest weight. 2 - the two vertices linked by this edge are merged: the new vertex is the parent of the two merged vertices 3 - the weight of the edges linking the new vertex to the remaining vertices of the graph are updated according to the user provided function (weight_function) 4 - repeat until a single edge remain
The initial weight of the edges (xedge_weights) and the callback (weight_function) determine the shape of the hierarchy.
The weight_function callback can be anything that defining the operator() and should follow the following pattern:
struct my_weighter { …
template<typename graph_t, typename neighbours_t> void operator()(const graph_t &g, // the current state of the graph index_t fusion_edge_index, // the edge between the two vertices being merged index_t new_region, // the new vertex in the graph index_t merged_region1, // the first vertex merged index_t merged_region2, // the second vertex merged neighbours_t &new_neighbours){ // list of edges to be weighted (see below) … for (auto &n: new_neighbours) { … n.new_edge_weight() = new_edge_value; // define the weight of this edge } }
Each element in the parameter new_neighbours represent an edge between the new vertex and another vertex of the graph. For each element of the list, the following methods are available:
neighbour_vertex(): the other vertex
num_edges(): returns 2 if both the two merged vertices add an edge linking themselves with neighbour_vertex() and 1 otherwise
first_edge_index(): the index of the edge linking one of the merged region to neighbour_vertex()
second_edge_index(): the index of the edge linking the other merged region to neighbour_vertex() (only if num_edges()==2)
new_edge_weight(): weight of the new edge (THIS HAS TO BE DEFINED IN THE WEIGHTING FUNCTION)
new_edge_index(): the index of the new edge: the weighting function will probably have to track new weight values
Example of weighting function: binary_partition_tree_min_linkage
- Template Parameters:
graph_t –
weighter –
T –
- Parameters:
graph –
xedge_weights –
weight_function –
- Returns:
a node weighted tree
-
template<typename
graph_t
, typenameT
>
autohg
::
binary_partition_tree_min_linkage
(const graph_t &graph, const xt::xexpression<T> &xedge_weights)¶ Binary partition tree, i.e.
the agglomerative clustering, with the minimum/single linkage rule.
Given a graph :math:
G=(V, E)
, with initial edge weights :math:w
, the distance :math:d(X,Y)
between any two clusters :math:X
and :math:Y
is.. math::
d(X,Y) = \min \{w(\{x,y\}) | x \in X, y \in Y, \{x,y\} \in E \}
Regions are then iteratively merged following the above distance (closest first) until a single region remains.
- Template Parameters:
graph_t –
T –
- Parameters:
graph –
xedge_weights –
- Returns:
a node weighted tree
-
template<typename
graph_t
, typenameT
>
autohg
::
binary_partition_tree_complete_linkage
(const graph_t &graph, const xt::xexpression<T> &xedge_weights)¶ Binary partition tree, i.e.
the agglomerative clustering, with the maximum/complete linkage rule.
Given a graph :math:
G=(V, E)
, with initial edge weights :math:w
, the distance :math:d(X,Y)
between any two clusters :math:X
and :math:Y
is.. math::
d(X,Y) = \max \{w(\{x,y\}) | x \in X, y \in Y, \{x,y\} \in E \}
Regions are then iteratively merged following the above distance (closest first) until a single region remains
- Template Parameters:
graph_t –
T –
- Parameters:
graph –
xedge_weights –
- Returns:
a node weighted tree
-
template<typename
graph_t
, typenameT
>
autohg
::
binary_partition_tree_average_linkage
(const graph_t &graph, const xt::xexpression<T> &xedge_weights, const xt::xexpression<T> &xedge_weight_weights)¶ \in E} w({x,y})`.
Regions are then iteratively merged following the above distance (closest first) until a single region remains
- Template Parameters:
graph_t –
T –
- Parameters:
graph –
xedge_weights –
xedge_weight_weights –
- Returns:
a node weighted tree
-
template<typename
graph_t
, typenameT
>
autohg
::
binary_partition_tree_exponential_linkage
(const graph_t &graph, const xt::xexpression<T> &xedge_weights, const typename T::value_type &alpha, const xt::xexpression<T> &xedge_weight_weights)¶ ))`.
Regions are then iteratively merged following the above distance (closest first) until a single region remains
Note that:
:math:
\\apha=0
is equivalent to average linkage clustering:math:
\\alpha=-\infty
is equivalent to single linkage clustering:math:
\\alpha=+\infty
is equivalent to complete linkage clustering
See:
Nishant Yadav, Ari Kobren, Nicholas Monath, Andrew Mccallum ; Supervised Hierarchical Clustering with Exponential Linkage Proceedings of the 36th International Conference on Machine Learning, PMLR 97:6973-6983, 2019.
- Template Parameters:
graph_t –
T –
- Parameters:
graph –
xedge_weights –
alpha –
xedge_weight_weights –
- Returns:
a node weighted tree
-
template<typename
graph_t
, typenameT1
, typenameT2
>
autohg
::
binary_partition_tree_ward_linkage
(const graph_t &graph, const xt::xexpression<T1> &xvertex_centroids, const xt::xexpression<T2> &xvertex_sizes, const std::string &altitude_correction = "max")¶ Binary partition tree, i.e.
the agglomerative clustering, with the Ward linkage rule.
Given a graph :math:
G=(V, E)
, with initial edge weights :math:w
with associated weights :math:w, the distance :math:
d(X,Y)between any two clusters :math:
Xand :math:
Y` is.. math::
d(X,Y) = \\frac{| X |\\times| Y |}{| X |+| Y |} \| \\vec{X} - \\vec{Y} \|^2
where :math:
\\vec{X}
and :math:\\vec{Y}
are the centroids of :math:X
and :math:Y
.Regions are then iteratively merged following the above distance (closest first) until a single region remains
Note that the Ward distance is not necessarily strictly increasing when processing a non complete graph. This can be corrected afterward with an altitude correction strategy. Valid values for
altitude correction
are:- ``"none"``: nothing is done and the altitude of a node is equal to the Ward distance between its 2 children; this may not be non-decreasing - ``"max"``: the altitude of a node :math:`n` is defined as the maximum of the the Ward distance associated to each node in the subtree rooted in :math:`n`.
- Template Parameters:
graph_t –
T1 –
T2 –
- Parameters:
graph –
xvertex_centroids – Centroids of the graph vertices (must be a 2d array)
xvertex_sizes – Size (number of elements) of the graph vertices
altitude_correction – can be
"none"
or"max"
(default)
- Returns:
a node weighted tree
-
template<typename
tree_t
, typenamealtitude_t
>
decltype(auto)hg
::
make_node_weighted_tree
(tree_t &&tree, altitude_t &&node_altitude)¶
-
template<typename
tree_t
, typenamenode_map_t
>
decltype(auto)hg
::
make_remapped_tree
(tree_t &&tree, node_map_t &&node_map)¶
-
template<typename
graph_t
, typenameT
>
autohg
::
component_tree_max_tree
(const graph_t &graph, const xt::xexpression<T> &xvertex_weights)¶ Construct the Max Tree of the vertex weighted graph.
The Min/Max Tree structure were proposed in [1], [2]. The algorithm used in this implementation was first described in [3].
[1] Ph. Salembier, A. Oliveras, and L. Garrido, “Anti-extensive connected operators for image
and sequence processing,” IEEE Trans. Image Process., vol. 7, no. 4, pp. 555-570, Apr. 1998.
[2] Ro. Jones, “Connected filtering and segmentation using component trees,” Comput. Vis. Image Understand., vol. 75, no. 3, pp. 215-228, Sep. 1999.
[3] Ch. Berger, T. Geraud, R. Levillain, N. Widynski, A. Baillard, and E. Bertin, “Effective
Component Tree Computation with Application to Pattern Recognition in Astronomical Imaging,” IEEE ICIP 2007.
- Template Parameters:
graph_t –
T –
- Parameters:
graph – input graph
vertex_weights – graph vertex weights
- Returns:
a node weighted tree
-
template<typename
graph_t
, typenameT
>
autohg
::
component_tree_min_tree
(const graph_t &graph, const xt::xexpression<T> &xvertex_weights)¶ Construct the Min Tree of the vertex weighted graph.
The Min/Max Tree structure were proposed in [1], [2]. The algorithm used in this implementation was first described in [3].
[1] Ph. Salembier, A. Oliveras, and L. Garrido, “Anti-extensive connected operators for image
and sequence processing,” IEEE Trans. Image Process., vol. 7, no. 4, pp. 555-570, Apr. 1998.
[2] Ro. Jones, “Connected filtering and segmentation using component trees,” Comput. Vis. Image Understand., vol. 75, no. 3, pp. 215-228, Sep. 1999.
[3] Ch. Berger, T. Geraud, R. Levillain, N. Widynski, A. Baillard, and E. Bertin, “Effective
Component Tree Computation with Application to Pattern Recognition in Astronomical Imaging,” IEEE ICIP 2007.
- Template Parameters:
graph_t –
T –
- Parameters:
graph – input graph
vertex_weights – graph vertex weights
- Returns:
a node weighted tree
-
template<typename
tree_t
, typenamealtitude_t
>
decltype(auto)hg
::
make_node_weighted_tree_and_mst
(tree_t &&tree, altitude_t &&node_altitude, array_1d<index_t> &&mst_edge_map)¶
-
template<typename
graph_t
, typenameT
>
autohg
::
bpt_canonical
(const graph_t &graph, const xt::xexpression<T> &xedge_weights)¶ Compute the canonical binary partition tree (or binary partition tree by altitude ordering) of the given edge weighted graph.
The algorithm returns a tuple composed of:
the binary partition tree,
the levels of the vertices of the tree,
the minimum spanning tree of the given graph that corresponds to this tree.
L. Najman, J. Cousty, B. Perret. Playing with Kruskal: algorithms for morphological trees in edge-weighted graphs. In, 11th International Symposium on Mathematical Morphology, ISMM 2013, Uppsala, Sweden, Mai 2013.
- Template Parameters:
graph_t –
T –
- Parameters:
graph –
xedge_weights –
- Returns:
-
template<typename
criterion_t
>
autohg
::
simplify_tree
(const tree &t, const criterion_t &criterion, bool process_leaves = false)¶ Creates a copy of the current Tree and deletes the nodes such that the criterion function is true.
Also returns an array that maps any node index i of the new tree, to the index of this node in the original tree.
The criterion function is a predicate that associates true (this node must be deleted) or false (do not delete this node) to a node index (with operator ()).
- Template Parameters:
criterion_t –
- Parameters:
t – input tree
criterion – For any vertex n of the tree, n has to be removed if criterion(n) == true
process_leaves – If false, a leaf vertex will never be removed disregarding the value of criterion.
- Returns:
-
template<typename
graph_t
, typenameT
>
autohg
::
quasi_flat_zone_hierarchy
(const graph_t &graph, const xt::xexpression<T> &xedge_weights)¶ Compute the quasi-flat zone hierarchy of an edge weighted graph.
For a given positive real value lamba:
a set of vertices X is lambda-connected if, for any two vertices x, y in X there exists an xy-path in X composed of edges of weights smaller of equal than lambda;
a lambda-connected component is a lambda-connected set of maximal extent;
the set of lambda-connected components forms a partition, called lambda-partition, of the graph vertices.
The quasi-flat zone hierarchy is composed of the sequence of lambda-partitions obtained for all lambda in edge_weights.
- Template Parameters:
graph_t – Input graph type
T – xepression derived type of input edge weights
- Parameters:
graph – Input graph
xedge_weights – Input graph weights
- Returns:
A node weighted tree
-
template<typename
graph_t
, typenametree_t
, typenameT
>
autohg
::
saliency_map
(const graph_t &graph, const tree_t &tree, const xt::xexpression<T> &xaltitudes)¶ Compute the saliency map of the given hierarchy for the given graph.
The saliency map is a weighting of the graph edges. The weight of an edge {x, y} is the altitude of the lowest common ancestor of x and y in the hierarchy.
- Template Parameters:
graph_t – Input graph type
tree_t – Input tree type
T – xepression derived type of input altitudes
- Parameters:
graph – Input graph
tree – Input tree
xaltitudes – Input node altitudes of the given tree
- Returns:
An array of shape (num_edges(graph)) and with the same value type as T.
-
template<typename
tree_t
>
autohg
::
tree_2_binary_tree
(const tree_t &tree)¶ Transforms a tree into a binary tree.
Each non-leaf node of the input tree must have at least 2 children!
Whenever a non-leaf node :math:
n
with :math:k > 2
children is found:an extra node :math:
m
is created;the first 2 children of :math:
n
become children of the new node :math:m
; andthe new node :math:
m
becomes the first child of :math:n
.
The number of children of :math:
n
is thus reduced by 1. This operation is repeated :math:k-2
times, i.e. until :math:n
has only 2 children.- Template Parameters:
tree_t – Input tree type
- Parameters:
tree – Input tree
- Returns:
-
template<typename
graph_t
, typenameT
, typenameF
>
autohg
::
watershed_hierarchy_by_attribute
(const graph_t &graph, const xt::xexpression<T> &xedge_weights, const F &attribute_functor)¶ Computes a hierarchical watershed for the given regional attribute.
The algorithm used is described in:
Laurent Najman, Jean Cousty, Benjamin Perret: Playing with Kruskal: Algorithms for Morphological Trees in Edge-Weighted Graphs. ISMM 2013: 135-146
The regional attribute is specified by an attribute functor, that is a function that
takes 2 input parameter: a binary partition tree, the altitudes of its nodes;
returns an 1d array giving the attribute value for each node of the input tree. The computed regional attribute must be scalar, positive and increasing (the attribute value of a node is smaller than or equal to the attribute value of its parent).
- Template Parameters:
graph_t –
T –
F –
- Parameters:
graph – input graph
xedge_weights – input graph edge weights
attribute_functor – function that computes the attribute value from a tree and its node altitudes
- Returns:
-
template<typename
graph_t
, typenameT1
, typenameT2
>
autohg
::
watershed_hierarchy_by_minima_ordering
(const graph_t &graph, const xt::xexpression<T1> &xedge_weights, const xt::xexpression<T2> &xminima_ranks)¶ Computes a hierarchical watershed for the given minima ordering.
The algorithm used is described in:
Laurent Najman, Jean Cousty, Benjamin Perret: Playing with Kruskal: Algorithms for Morphological Trees in Edge-Weighted Graphs. ISMM 2013: 135-146
The ranking ranking of the minima of the given edge weighted graph (G,w) is given as vertex weights with values in {0..n} with n the number of minima of (G,w). It must satisfy the following pre-conditions:
each minimum of (G,w) contains at least one non zero vertex,
all non zero vertices in a minimum have the same weight,
there is no non zero value vertex outside minima, and
no two minima contain non zero vertices with the same weight.
The altitude associated to each minimum is a non decreasing 1d array of size n + 1 with non negative values. Note that the first entry of the minima altitudes array, ie. the value at index 0, does not represent a minimum and its value should be 0.
- Template Parameters:
graph_t –
T1 –
T2 –
- Parameters:
graph – input graph
xedge_weights – input graph edge weights
xminima_ranks – input graph vertex weights containing the rank of each minima of the input edge weighted graph
- Returns:
-
template<typename
graph_t
, typenameT1
, typenameT2
>
autohg
::
watershed_hierarchy_by_area
(const graph_t &graph, const xt::xexpression<T1> &edge_weights, const xt::xexpression<T2> &vertex_area)¶
-
template<typename
graph_t
, typenameT1
>
autohg
::
watershed_hierarchy_by_area
(const graph_t &graph, const xt::xexpression<T1> &xedge_weights)¶
-
template<typename
graph_t
, typenameT1
, typenameT2
>
autohg
::
watershed_hierarchy_by_volume
(const graph_t &graph, const xt::xexpression<T1> &edge_weights, const xt::xexpression<T2> &vertex_area)¶
-
template<typename
graph_t
, typenameT1
>
autohg
::
watershed_hierarchy_by_volume
(const graph_t &graph, const xt::xexpression<T1> &xedge_weights)¶
-
template<typename
graph_t
, typenameT
>
autohg
::
watershed_hierarchy_by_dynamics
(const graph_t &graph, const xt::xexpression<T> &edge_weights)¶
-
template<typename
graph_t
, typenameT
>
autohg
::
fit_contour_2d
(const graph_t &graph, const embedding_grid_2d &embedding, const xt::xexpression<T> &xedge_weights)¶ Construct a contour_2d object from a graph cut of a 2d image with a 4 adjacency (non zero edges are part of the cut).
- Template Parameters:
graph_t –
T –
- Parameters:
graph –
embedding –
xedge_weights –
- Returns:
-
template<typename
rag_graph_t
, typenamegraph_t
, typenameT
>
autohg
::
rag_2d_vertex_perimeter_and_edge_length
(const rag_graph_t &rag_graph, const xt::xexpression<T> &xvertex_map, const xt::xexpression<T> &xedge_map, const embedding_grid_2d &embedding, const graph_t &graph, double epsilon = 0.1, bool relative_epsilon = true, int min_size = 2)¶ Estimate the vertex perimeter and the length of the frontier associated to the edges of a region adjacency graph constructed on a 2d 4 adjacency graph.
The region boundaries are simplified with Ramer–Douglas–Peucker algorithm and is controlled by the parameters epsilon, relative_epsilon, min_size. See function subdivide of the class contour_2d for more information.
- Template Parameters:
graph_t –
- Parameters:
rag_graph – Region Adjacency Graph
xvertex_map – Vertex map of the rag_graph
xedge_map – Edge map of the rag_graph
embedding – 2d shape of the input graph
graph – input graph on which the region adjacency graph has been build: must be a 4 adjacency graph whose shape correspond to the given embedding
epsilon – larger epsilon values will provide stronger contour shapes simplification
relative_epsilon – Is epsilon given in relative or absolute units
min_size – Boundaries elements smaller than min_size will be deleted
- Returns:
a pair composed of two 1d arrays: vertex_perimeter and edge_length.
-
template<typename
graph_t
>
autohg
::
rag_2d_vertex_perimeter_and_edge_length
(const region_adjacency_graph &rag, const embedding_grid_2d &embedding, const graph_t &graph, double epsilon = 0.1, bool relative_epsilon = true, int min_size = 2)¶ Estimate the vertex perimeter and the length of the frontier associated to the edges of a region adjacency graph constructed on a 2d 4 adjacency graph.
The region boundaries are simplified with Ramer–Douglas–Peucker algorithm and is controlled by the parameters epsilon, relative_epsilon, min_size. See function subdivide of the class contour_2d for more information.
- Template Parameters:
graph_t –
- Parameters:
rag – Region Adjacency Graph
embedding – 2d shape of the input graph
graph – input graph on which the region adjacency graph has been build: must be a 4 adjacency graph whose shape correspond to the given embedding
epsilon – larger epsilon values will provide stronger contour shapes simplification
relative_epsilon – Is epsilon given in relative or absolute units
min_size – Boundaries elements smaller than min_size will be deleted
- Returns:
a pair composed of two 1d arrays: vertex_perimeter and edge_length.
-
inline auto
hg
::
get_4_adjacency_implicit_graph
(const embedding_grid_2d &embedding)¶ Create a 4 adjacency implicit regular graph for the given embedding.
- Parameters:
embedding –
- Returns:
-
inline auto
hg
::
get_8_adjacency_implicit_graph
(const embedding_grid_2d &embedding)¶ Create of 4 adjacency implicit regular graph for the given embedding.
- Parameters:
embedding –
- Returns:
-
inline auto
hg
::
get_4_adjacency_graph
(const embedding_grid_2d &embedding)¶ Create of 4 adjacency explicit regular graph for the given embedding.
- Parameters:
embedding –
- Returns:
-
inline auto
hg
::
get_8_adjacency_graph
(const embedding_grid_2d &embedding)¶ Create of 8 adjacency explicit regular graph for the given embedding.
- Parameters:
embedding –
- Returns:
-
template<typename
graph_t
, typenameT
, typenameresult_type
= typename T::value_type>
autohg
::
graph_4_adjacency_2_khalimsky
(const graph_t &graph, const embedding_grid_2d &embedding, const xt::xexpression<T> &xedge_weights, bool add_extra_border = false, result_type extra_border_value = 0)¶ Represents a 4 adjacency edge weighted regular graph in 2d Khalimsky space.
- Parameters:
embedding –
- Returns:
-
template<typename
T
, typenameresult_type
= typename T::value_type>
autohg
::
khalimsky_2_graph_4_adjacency
(const xt::xexpression<T> &xkhalimsky, bool extra_border = false)¶ Transforms a contour map represented in 2d Khalimsky space into a weighted 4 adjacency edge weighted regular graph (0-face and 2-face of the Khalimsky space are ignored).
- Parameters:
embedding –
- Returns:
-
template<typename
graph_t
, typenameT1
, typenameT2
>
autohg
::
oriented_watershed
(const graph_t &graph, const embedding_grid_2d &embedding, const xt::xexpression<T1> &xedge_weights, const xt::xexpression<T2> &xedge_orientations = array_nd<int>())¶ Compute the oriented watershed as described in [ArbelaezPAMI2011]_ .
Given a 4 adjacency graph with edge boundary probabilities and estimated boundary orientations, the algorithms computes:
a region adjacency graph of the watershed regions of the edge boundary probabilities
the boundaries between watershed regions are vectorized and simplified (see contour_2d class)
the orientation of each boundary element is estimated
the edge boundary probabilities are reweighted according to the concordance between user provided boundary orientations and estimated orientation of boundary elements
the weight of the region adjacency graph edges as the mean value of reweighted edge boundary probabilities on the frontier between the 2 regions
The algorithm returns the region adjacency graph of watershed pixels and its edge weights.
.. [ArbelaezPAMI2011] Arbelaez, P., Maire, M., Fowlkes, C., & Malik, J.. Contour detection and hierarchical image segmentation. IEEE transactions on pattern analysis and machine intelligence, 33(5), 898-916.
- Template Parameters:
graph_t –
T1 –
T2 –
- Parameters:
graph –
embedding –
xedge_weights –
xedge_orientations –
- Returns:
-
template<typename
graph_t
, typenameT1
, typenameT2
>
autohg
::
mean_pb_hierarchy
(const graph_t &graph, const embedding_grid_2d &embedding, const xt::xexpression<T1> &xedge_weights, const xt::xexpression<T2> &xedge_orientations = array_nd<int>())¶ Compute the mean probability boundary hierarchy as described in [ArbelaezPAMI2011]_ .
Given a 4 adjacency graph with edge boundary probabilities and estimated boundary orientations, the algorithms computes:
the oriented watershed of the given graph
the average linkage clustering ot the oriented watershed
The algorithm returns the region adjacency graph of watershed pixels and teh valued tree computed on this graph.
.. [ArbelaezPAMI2011] Arbelaez, P., Maire, M., Fowlkes, C., & Malik, J.. Contour detection and hierarchical image segmentation. IEEE transactions on pattern analysis and machine intelligence, 33(5), 898-916.
- Template Parameters:
graph_t –
T1 –
T2 –
- Parameters:
graph –
embedding –
xedge_weights –
xedge_orientations –
- Returns:
-
template<typename
T
>
autohg
::
component_tree_tree_of_shapes_image2d
(const xt::xexpression<T> &ximage, tos_padding padding = tos_padding::mean, bool original_size = true, bool immersion = true, index_t exterior_vertex = 0)¶ Computes the tree of shapes of a 2d image.
The Tree of Shapes was described in [1].
The algorithm used in this implementation was first described in [2].
The tree is computed in the interpolated multivalued Khalimsky space to provide a continuous and autodual representation of input image.
If padding is different from tos_padding::none, an extra border of pixels is added to the input image before anything else. This will ensure the existence of a shape encompassing all the shapes inside the input image (if exterior_vertex is inside the extra border): this shape will be the root of the tree. The padding value can be:
0 is padding == tos_padding::zero
the mean value of the boundary pixels of the input image if padding == tos_padding::mean
If original_size is true, all the nodes corresponding to pixels not belonging to the input image are removed (except for the root node). If original_size is false, the returned tree is the tree constructed in the interpolated/padded space. In practice if the size of the input image is (h, w), the leaves of the returned tree will correspond to an image of size:
(h, w) if original_size is true;
(h * 2 - 1, w * 2 - 1) is original_size is false and padding is tos_padding::none; and
((h + 2) * 2 - 1, (w + 2) * 2 - 1) otherwise.
:Advanced options:
Use with care the following options may lead to unexpected results:
Immersion defines if the initial image should be first converted as an equivalent continuous representation called a plain map. If the immersion is deactivated the level lines of the shapes of the image may intersect (if the image is not well composed) and the result of the algorithm is undefined. If immersion is deactivated, the factor :math:
*2 - 1
has to be removed in the result sizes given above.Exterior_vertex defines the linear coordinates of the pixel corresponding to the exterior (interior and exterior of a shape is defined with respect to this point). The coordinate of this point must be given in the padded/interpolated space.
[1] Pa. Monasse, and F. Guichard, “Fast computation of a contrast-invariant image representation,” Image Processing, IEEE Transactions on, vol.9, no.5, pp.860-872, May 2000
[2] Th. Géraud, E. Carlinet, S. Crozet, and L. Najman, “A Quasi-linear Algorithm to Compute the Tree
of Shapes of nD Images”, ISMM 2013.
- Template Parameters:
T –
- Parameters:
ximage – Must be a 2d array
padding – Defines if an extra boundary of pixels is added to the original image (see enum tos_padding).
original_size – remove all nodes corresponding to interpolated/padded pixels
exterior_vertex – linear coordinate of the exterior point
- Returns:
a node weighted tree
-
inline auto
hg
::
read_pink_graph
(std::istream &in)¶
-
inline auto
hg
::
read_pink_graph
(const std::string &filename)¶
-
template<typename
graph_t
, typenameT1
, typenameT2
, typenameS
>
autohg
::
save_pink_graph
(std::ostream &out, const graph_t &graph, const xt::xexpression<T1> &xvertex_values = xt::xscalar<char>(0), const xt::xexpression<T2> &xedge_values = xt::xscalar<char>(0), S &shape = std::vector<std::size_t>())¶
-
template<typename
graph_t
, typenameT1
, typenameT2
, typenameS
>
voidhg
::
save_pink_graph
(const std::string &filename, const graph_t &graph, const xt::xexpression<T1> &xvertex_values = xt::xscalar<char>(0), const xt::xexpression<T2> &xedge_values = xt::xscalar<char>(0), S &shape = std::vector<std::size_t>())¶
-
inline array_nd<unsigned char>
hg
::
read_image_pnm
(const char *filename)¶ Read the given pnm image (pbm, pgm or ppm formats).
Current the following pnm specification are supported P1 binary ascii: supported P2 byte ascii: supported (max value <= 255) P3 RGB ascii: supported (max value <= 255) P4 binary raw: NOT supported P5 byte raw: supported (max value <= 255) P6 RGB rax: supported (max value <= 255)
- Parameters:
filename – Path to the file to read
- Returns:
an array with pixel data
-
template<typename
T
>
voidhg
::
save_image_pnm
(const char *filename, const xt::xexpression<T> &ximage)¶ Save an array as a pnm file (pgm or ppm).
The array value_type MUST be unsigned char. If the array has 2 dimensions it is saved as a pgm raw file (format P5). If the array has 3 dimensions, the size of the third dimension must be 3 and it is saved as a ppm raw file (format P6).
If the provided filename already exists, it will be overwritten !
- Template Parameters:
T – xexpression derived type for ximage
- Parameters:
filename – Path to the file to write.
ximage – The array containing the data to save.
-
inline auto
hg
::
read_tree
(std::istream &in)¶
-
template<typename
RandomAccessIterator
, typenameCompare
>
voidhg
::
stable_sort
(RandomAccessIterator xs, RandomAccessIterator xe, Compare comp)¶
-
template<typename
RandomAccessIterator
, typenameCompare
>
voidhg
::
sort
(RandomAccessIterator xs, RandomAccessIterator xe, Compare comp)¶
-
template<typename
RandomAccessIterator
>
voidhg
::
stable_sort
(RandomAccessIterator xs, RandomAccessIterator xe)¶
-
template<typename
T
, typenameCompare
>
voidhg
::
stable_sort
(xt::xexpression<T> &arrayx, Compare comp)¶
-
template<typename
RandomAccessIterator
>
voidhg
::
sort
(RandomAccessIterator xs, RandomAccessIterator xe)¶
-
template<typename
T
, typenameCompare
>
autohg
::
arg_sort
(const xt::xexpression<T> &arrayx, Compare comp)¶
-
template<typename
T
, typenameCompare
>
autohg
::
stable_arg_sort
(const xt::xexpression<T> &arrayx, Compare comp)¶
-
template<typename
graph_t
>
auto &hg
::
source
(const indexed_edge<typename graph::graph_traits<graph_t>::vertex_descriptor, typename graph::graph_traits<graph_t>::edge_index> &e, const graph_t&)¶ Source vertex of an edge.
- Template Parameters:
graph_t –
- Parameters:
e –
- Returns:
-
template<typename
graph_t
>
auto &hg
::
target
(const indexed_edge<typename graph::graph_traits<graph_t>::vertex_descriptor, typename graph::graph_traits<graph_t>::edge_index> &e, const graph_t&)¶ Target vertex of an edge.
- Template Parameters:
graph_t –
- Parameters:
e –
- Returns:
-
template<typename
graph_t
>
auto &hg
::
index
(const indexed_edge<typename graph::graph_traits<graph_t>::vertex_descriptor, typename graph::graph_traits<graph_t>::edge_index> &e, const graph_t&)¶ Index of an edge.
- Template Parameters:
graph_t –
- Parameters:
e –
- Returns:
-
template<bool
vectorial
= true, typenameT
>
autohg
::
make_light_axis_view
(T &e, index_t position = 0)¶
-
template<typename
embedding_t
>
std::pair<typename hg::regular_graph<embedding_t>::out_edge_iterator, typename hg::regular_graph<embedding_t>::out_edge_iterator>hg
::
out_edges
(typename hg::regular_graph<embedding_t>::vertex_descriptor u, const hg::regular_graph<embedding_t> &g)¶
-
template<typename
embedding_t
>
std::pair<typename hg::regular_graph<embedding_t>::out_edge_iterator, typename hg::regular_graph<embedding_t>::out_edge_iterator>hg
::
in_edges
(typename hg::regular_graph<embedding_t>::vertex_descriptor u, const hg::regular_graph<embedding_t> &g)¶
-
template<typename
embedding_t
>
hg::regular_graph<embedding_t>::vertices_size_typehg
::
num_vertices
(const hg::regular_graph<embedding_t> &g)¶
-
template<typename
embedding_t
>
std::pair<typename hg::regular_graph<embedding_t>::vertex_iterator, typename hg::regular_graph<embedding_t>::vertex_iterator>hg
::
vertices
(const hg::regular_graph<embedding_t> &g)¶
-
template<typename
embedding_t
>
std::pair<typename hg::regular_graph<embedding_t>::adjacency_iterator, typename hg::regular_graph<embedding_t>::adjacency_iterator>hg
::
adjacent_vertices
(typename hg::regular_graph<embedding_t>::vertex_descriptor u, const hg::regular_graph<embedding_t> &g)¶
-
template<typename
embedding_t
>
hg::regular_graph<embedding_t>::degree_size_typehg
::
out_degree
(const typename hg::regular_graph<embedding_t>::vertex_descriptor v, const hg::regular_graph<embedding_t> &g)¶
-
template<typename
embedding_t
>
hg::regular_graph<embedding_t>::degree_size_typehg
::
in_degree
(const typename hg::regular_graph<embedding_t>::vertex_descriptor v, const hg::regular_graph<embedding_t> &g)¶
-
template<typename
embedding_t
>
hg::regular_graph<embedding_t>::degree_size_typehg
::
degree
(const typename hg::regular_graph<embedding_t>::vertex_descriptor v, const hg::regular_graph<embedding_t> &g)¶
-
inline auto
hg
::
leaves_to_root_iterator
(const tree &t, leaves_it leaves_opt = leaves_it::include, root_it root_opt = root_it::include)¶
-
inline auto
hg
::
root_to_leaves_iterator
(const tree &t, leaves_it leaves_opt = leaves_it::include, root_it root_opt = root_it::include)¶
-
inline std::pair<tree::children_iterator, tree::children_iterator>
hg
::
children
(const tree::vertex_descriptor v, const tree &g)¶
-
inline hg::tree::degree_size_type
hg
::
degree
(const typename hg::tree::vertex_descriptor v, const hg::tree &g)¶
-
inline hg::tree::degree_size_type
hg
::
in_degree
(const typename hg::tree::vertex_descriptor v, const hg::tree &g)¶
-
inline hg::tree::degree_size_type
hg
::
out_degree
(const typename hg::tree::vertex_descriptor v, const hg::tree &g)¶
-
inline std::pair<typename hg::tree::vertex_iterator, typename hg::tree::vertex_iterator>
hg
::
vertices
(const hg::tree &g)¶
-
inline std::pair<typename hg::tree::edge_iterator, typename hg::tree::edge_iterator>
hg
::
edges
(const hg::tree &g)¶
-
inline std::pair<typename hg::tree::adjacency_iterator, typename hg::tree::adjacency_iterator>
hg
::
adjacent_vertices
(typename hg::tree::vertex_descriptor v, const hg::tree &g)¶
-
inline std::pair<hg::tree::out_edge_iterator, hg::tree::out_edge_iterator>
hg
::
out_edges
(hg::tree::vertex_descriptor v, const hg::tree &g)¶
-
inline std::pair<hg::tree::out_edge_iterator, hg::tree::out_edge_iterator>
hg
::
in_edges
(hg::tree::vertex_descriptor v, const hg::tree &g)¶
-
template<typename
T
>
autohg
::
find_region
(tree::vertex_descriptor v, const typename T::value_type &lambda, const T &altitudes, const tree &tree)¶
-
template<typename
T1
, typenameT2
, typenameT3
>
autohg
::
find_region
(const xt::xexpression<T1> &xvertices, const xt::xexpression<T2> &xlambdas, const xt::xexpression<T3> &xaltitudes, const tree &t)¶
-
inline auto
hg
::
lowest_common_ancestor
(tree::vertex_descriptor v1, tree::vertex_descriptor v2, const tree &t)¶
-
template<typename
T
>
autohg
::
lowest_common_ancestor
(const xt::xexpression<T> &xvertices_1, const xt::xexpression<T> &xvertices_2, const tree &t)¶
-
template<typename
T
>
const auto &hg
::
edge_from_index
(const typename undirected_graph<T>::vertex_descriptor v, const undirected_graph<T> &g)¶
-
template<typename
T
>
hg::undirected_graph<T>::vertices_size_typehg
::
num_vertices
(const hg::undirected_graph<T> &g)¶
-
template<typename
T
>
hg::undirected_graph<T>::edges_size_typehg
::
num_edges
(const hg::undirected_graph<T> &g)¶
-
template<typename
T
>
hg::undirected_graph<T>::degree_size_typehg
::
degree
(typename hg::undirected_graph<T>::vertex_descriptor v, const hg::undirected_graph<T> &g)¶
-
template<typename
T
>
hg::undirected_graph<T>::degree_size_typehg
::
in_degree
(typename hg::undirected_graph<T>::vertex_descriptor v, const hg::undirected_graph<T> &g)¶
-
template<typename
T
>
hg::undirected_graph<T>::degree_size_typehg
::
out_degree
(typename hg::undirected_graph<T>::vertex_descriptor v, const hg::undirected_graph<T> &g)¶
-
template<typename
T
>
hg::undirected_graph<T>::vertex_descriptorhg
::
add_vertex
(hg::undirected_graph<T> &g)¶
-
template<typename
T
>
voidhg
::
add_vertices
(size_t num, hg::undirected_graph<T> &g)¶
-
template<typename
T
>
hg::undirected_graph<T>::edge_descriptorhg
::
add_edge
(typename hg::undirected_graph<T>::vertex_descriptor v1, typename hg::undirected_graph<T>::vertex_descriptor v2, hg::undirected_graph<T> &g)¶
-
template<typename
T
>
voidhg
::
remove_edge
(typename hg::undirected_graph<T>::edge_index_t ei, hg::undirected_graph<T> &g)¶
-
template<typename
T
>
voidhg
::
set_edge
(typename hg::undirected_graph<T>::edge_index_t ei, typename hg::undirected_graph<T>::vertex_descriptor v1, typename hg::undirected_graph<T>::vertex_descriptor v2, hg::undirected_graph<T> &g)¶
-
template<typename
T
>
std::pair<typename hg::undirected_graph<T>::vertex_iterator, typename hg::undirected_graph<T>::vertex_iterator>hg
::
vertices
(const hg::undirected_graph<T> &g)¶
-
template<typename
T
>
std::pair<typename hg::undirected_graph<T>::edge_iterator, typename hg::undirected_graph<T>::edge_iterator>hg
::
edges
(const hg::undirected_graph<T> &g)¶
-
template<typename
T
>
std::pair<typename hg::undirected_graph<T>::out_edge_iterator, typename hg::undirected_graph<T>::out_edge_iterator>hg
::
out_edges
(typename hg::undirected_graph<T>::vertex_descriptor v, const hg::undirected_graph<T> &g)¶
-
template<typename
T
>
std::pair<typename hg::undirected_graph<T>::out_edge_iterator, typename hg::undirected_graph<T>::out_edge_iterator>hg
::
in_edges
(typename hg::undirected_graph<T>::vertex_descriptor v, const hg::undirected_graph<T> &g)¶
-
template<typename
T
>
std::pair<typename hg::undirected_graph<T>::adjacency_iterator, typename hg::undirected_graph<T>::adjacency_iterator>hg
::
adjacent_vertices
(typename hg::undirected_graph<T>::vertex_descriptor v, const hg::undirected_graph<T> &g)¶
-
template<typename
lambda_t
>
voidhg
::
parfor
(index_t start_index, index_t end_index, lambda_t fun, index_t step_size = 1)¶
-
template<typename
T1
, typenameT2
>
voidhg
::
extend
(T1 &a, const T2 &b)¶ Insert all elements of collection b at the end of collection a.
- Template Parameters:
T1 – must have an insert method (STL like) and a range interface (begin, end)
T2 – must have a range interface (begin, end)
- Parameters:
a –
b –
-
inline void
hg
::
unreachable
()¶
-
struct
hg::hg
::
accumulator_argmax
¶ Public Functions
-
struct
hg::hg
::
accumulator_argmin
¶ Public Functions
-
struct
hg::hg
::
accumulator_counter
¶ Public Functions
-
struct
hg::hg
::
accumulator_first
¶ Public Functions
-
struct
hg::hg
::
accumulator_last
¶ Public Functions
-
struct
hg::hg
::
accumulator_max
¶ Public Functions
Public Static Functions
-
template<typename
value_type
>
static inline value_typeinit_value
()¶
-
template<typename
value_type
>
static inline value_typereduce
(const value_type &v1, const value_type &v2)¶
-
template<typename
-
struct
hg::hg
::
accumulator_mean
¶ Public Functions
-
struct
hg::hg
::
accumulator_min
¶ Public Functions
Public Static Functions
-
template<typename
value_type
>
static inline value_typeinit_value
()¶
-
template<typename
value_type
>
static inline value_typereduce
(const value_type &v1, const value_type &v2)¶
-
template<typename
-
struct
hg::hg
::
accumulator_prod
¶ Public Functions
Public Static Functions
-
template<typename
value_type
>
static inline value_typeinit_value
()¶
-
template<typename
value_type
>
static inline value_typereduce
(const value_type &v1, const value_type &v2)¶
-
template<typename
-
struct
hg::hg
::
accumulator_sum
¶ Public Functions
Public Static Functions
-
template<typename
value_type
>
static inline constexpr value_typeinit_value
()¶
-
template<typename
value_type
>
static inline value_typereduce
(const value_type &v1, const value_type &v2)¶
-
template<typename
-
class
hg::hg
::
assesser_fragmentation_optimal_cut
¶ - #include <fragmentation_curve.hpp>
This class is used to assess the optimal cuts of a hierarchy of partitions with respect to a given ground-truth labelisation of its base graph and the BCE measure.
Public Functions
-
template<typename
tree_t
, typenameT
>
inlineassesser_fragmentation_optimal_cut
(const tree_t &tree, const xt::xexpression<T> &xground_truth, optimal_cut_measure measure, const array_1d<index_t> &vertex_map = {}, size_t max_regions = 200)¶ Create an assesser for hierarchy optimal cuts w.r.t.
a given ground-truth partition of hierarchy leaves and the BCE quality measure. The algorithms will explore optimal cuts containing at most max_regions regions.
The ground truth labelisation must be normalized (i.e. its labels must be positive integers in the interval [0, num_regions[).
- Template Parameters:
tree_t – tree type
T – type of labels
- Parameters:
tree – input hierarchy
xground_truth – ground truth labelisation of the tree leaves
vertex_map – super-vertices map (if tree is built on a rag, leave empty otherwise)
max_regions – maximum number of regions in the considered cuts.
-
inline auto
fragmentation_curve
() const¶ Fragmentation curve, i.e.
for each number of region k between 1 and max_regions, the BCE score of the optimal cut with k regions.
- Returns:
-
inline auto
optimal_number_of_regions
() const¶ Number of regions in the optimal cut.
- Returns:
-
inline auto
optimal_score
() const¶ Score of the optimal cut.
- Returns:
-
inline auto
optimal_partition
(size_t num_regions = 0) const¶ Labelisation of the base graph that corresponds to the optimal cut with the given number of regions.
If the number of regions is equal to 0 (default), the global optimal cut it returned (it will contain get_optimal_number_of_regions regions).
- Parameters:
num_regions –
- Returns:
-
inline auto
straightened_altitudes
(bool gain_only = false, bool normalize_result = true) const¶ Compute tree node altitudes such that the horizontal cut of the resulting vertex valued hierarchy corresponds to the optimal cut of the tree.
- Parameters:
gain_only – If true, optimal cuts with a number of regions greater than optimal_number_of_regions() are ignored
normalize_result – Ensure that altitudes[tree(root)] == optimal_score
- Returns:
-
template<typename
-
template<typename
T
>
structCOMPILE_ERROR
¶ - #include <utils.hpp>
Do not use except if you want a compile error showing the type of the provided template parameter !
- Template Parameters:
T –
-
template<typename
value_type
= index_t>
structhg::hg
::
counting_iterator
: public hg::random_iterator_facade<counting_iterator<index_t>, index_t>¶ - #include <iterators.hpp>
Quick implementation of a counting iterator.
Counting iterator has an integer value and a step size. Incrementing the iterator simply add the step size to the current value of the iterator.
- Template Parameters:
value_type –
Public Types
-
using
self_type
= counting_iterator<value_type>¶
Public Functions
-
inline
counting_iterator
()¶
-
inline
counting_iterator
(value_type position, value_type step = 1)¶
-
inline void
increment
()¶
-
inline void
decrement
()¶
-
inline bool
equal
(const counting_iterator<value_type> &other) const¶
-
inline value_type
dereference
() const¶
-
template<typename
derived_type
, typenamevalue_t
, typenamereference_t
= value_t>
structhg::hg
::
forward_iterator_facade
¶ - #include <iterators.hpp>
Facade to ease the declaration of new forward iterators (quick reimplementation of boost facade iterator).
The iterator must inherit of the facade and declare the following 3 functions:
deference: obtain the element at the current position of the iterator
increment: move the iterator to next element
equal: test if two iterators are equal
- Template Parameters:
derived_type –
value_t –
reference_t –
Public Types
-
using
self_type
= forward_iterator_facade<derived_type, value_t, reference_t>¶
-
using
difference_type
= std::ptrdiff_t¶
-
using
reference
= reference_t¶
-
using
iterator_category
= std::forward_iterator_tag¶
Friends
-
inline friend reference
operator*
(const derived_type &lhs)¶
-
inline friend const friend reference operator* (derived_type &lhs)
-
inline friend derived_type
operator++
(derived_type &lhs)¶
-
inline friend derived_type
operator++
(derived_type &lhs, int)¶
-
inline friend bool
operator==
(const derived_type &lhs, const derived_type &rhs)¶
-
inline friend bool
operator!=
(const derived_type &lhs, const derived_type &rhs)¶
-
template<typename
value_t
= double>
structhg::hg
::
fragmentation_curve
¶ Public Functions
-
inline
fragmentation_curve
(array_1d<value_t> &&num_regions, array_1d<value_t> &&scores, size_t num_regions_ground_truth)¶
-
inline auto
num_regions_ground_truth
() const¶ Number of regions in the ground truth labelisation of the hierarchy base graph.
- Returns:
-
inline auto
optimal_number_of_regions
() const¶ Number of regions in the optimal cut.
- Returns:
-
inline auto
optimal_score
() const¶ Score of the optimal cut.
- Returns:
-
inline const auto &
scores
() const¶
-
inline const auto &
num_regions
() const¶
-
inline auto
num_regions_normalized
() const¶
-
inline
-
class
hg::hg
::
hierarchy_aligner
¶ - #include <alignment.hpp>
This class allows to project hierarchies build from coarse supervertices onto fine supervertices.
The class is contructed by providing a fine supervertices decomposition of a graph. Then the functions align_hierarchy allows to project a hierarchy, given as a tree or as a saliency map, onto the fine supervertices.
Given:
a graph g
a fine labelisation l1 of the vertices of g;
a tree t on g whose supervertices corresponds to the coarse labelisation l2 of the vertices of g; and
the altitudes a of the nodes of t. Let us denote:
given a vertex x of g and a labelisation l, l(x) is the region of l that contains x
given a region r of l1, s(r, l2) is the region R of l2 that has the largest intersection with r: s(r, l2) = arg_max(R in l2) |R \cap r| The projection of t onto l1 is a hierarchy given by the saliency map sm on g defined by: for all {x,y} in edges(g), sm({x,y}) = a(lca_t(s(l1(x), l2), s(l1(y), l2)))
See the following helper functions for instanciation
make_hierarchy_aligner_from_graph_cut
make_hierarchy_aligner_from_labelisation
make_hierarchy_aligner_from_hierarchy
Public Functions
-
inline
hierarchy_aligner
(region_adjacency_graph &&rag)¶
-
template<typename
T
>
inline autoalign_hierarchy
(const hg::tree &tree, const xt::xexpression<T> &xaltitudes) const¶
-
template<typename
tree_t
, typenamevalue_t
>
classhg::hg
::
horizontal_cut_explorer
¶
-
template<typename
value_t
>
structhg::hg
::
horizontal_cut_nodes
¶
-
template<typename
vertex_descriptor
, typenameedge_index_t
>
classhg::hg
::
indexed_edge
¶ - #include <indexed_edge.hpp>
An edge with a source vertex, a target vertex, and an index.
- Template Parameters:
vertex_descriptor –
edge_index_t –
Public Functions
-
inline
indexed_edge
(vertex_descriptor _source, vertex_descriptor _target, edge_index_t _index)¶
-
inline
indexed_edge
(const std::pair<vertex_descriptor, vertex_descriptor> &edge, edge_index_t _index)¶
-
inline
operator edge_index_t
() const¶
-
inline
operator std::pair<vertex_descriptor, vertex_descriptor>
() const¶
Public Members
-
vertex_descriptor
source
¶
-
vertex_descriptor
first
¶
-
union hg::indexed_edge::[anonymous] [anonymous]¶
-
vertex_descriptor
target
¶
-
vertex_descriptor
second
¶
-
union hg::indexed_edge::[anonymous] [anonymous]¶
-
edge_index_t
index
¶
-
template<typename
value_type
= index_t>
structhg::hg
::
irange
¶ - #include <iterators.hpp>
A range of integer value with a start value, a step size, and a maximum value (excluded).
irange(start, stop, step) is the equivalent to
python range(start, stop, step)
matlab start:step:stop
- Template Parameters:
value_type –
Public Functions
-
inline
irange
()¶
-
inline
irange
(value_type start, value_type stop, value_type step = 1)¶
-
inline auto
begin
() const¶
-
inline auto
end
() const¶
-
template<typename
iterator_t
>
structhg::hg
::
iterator_wrapper
¶ - #include <graph.hpp>
Simple wrapper over two iterators to create a “range” object usable in foreach loops.
- Template Parameters:
iterator_t –
Public Functions
-
inline
iterator_wrapper
(iterator_t &_first, iterator_t &_last)¶
-
inline
iterator_wrapper
(const std::pair<iterator_t, iterator_t> &p)¶
-
inline iterator_t
begin
()¶
-
inline iterator_t
end
()¶
-
inline iterator_t
begin
() const¶
-
inline iterator_t
end
() const¶
-
struct
hg::hg
::
logger
¶ Public Types
-
using
callback_list
= std::vector<std::function<void(const std::string&)>>¶
Public Static Functions
-
static inline bool &
trace_enabled
()¶
-
static inline callback_list &
callbacks
()¶
Public Static Attributes
-
static const std::size_t
MAX_MSG_SIZE
= 8096¶
-
using
-
template<typename
mst_t
>
structhg::hg
::
minimum_spanning_tree_result
¶ - #include <graph_core.hpp>
A simple structure to hold the result of minimum_spanning_tree function.
The structures holds 2 elements:
the minimum spanning tree (mst)
a map (mst_edge_map) that indicates for each edge of the mst, the corresponding edge index in the original graph
- Template Parameters:
mst_t –
-
template<typename
tree_t
, typenamealtitude_t
>
structhg::hg
::
node_weighted_tree
¶ - #include <common.hpp>
A simple structure to hold the result of hierarchy construction algorithms.
See make_node_weighted_tree for construction
- Template Parameters:
tree_t –
altitude_t –
-
template<typename
tree_t
, typenamealtitude_t
>
structhg::hg
::
node_weighted_tree_and_mst
¶ - #include <hierarchy_core.hpp>
A simple structure to hold the result of canonical bpt function.
See make_node_weighted_tree_and_mst for construction
- Template Parameters:
tree_t –
altitude_t –
mst_t –
-
template<typename
derived_type
, typenamevalue_t
, typenamereference_t
= value_t>
structhg::hg
::
random_iterator_facade
¶ - #include <iterators.hpp>
Facade to ease the declaration of new random iterators (quick reimplementation of boost facade iterator).
The iterator must inherit of the facade and declare the following 3 functions:
deference: obtain the element at the current position of the iterator
increment: move the iterator to next element
decrement: move the iterator to previous element
advance(n): Advance by n positions
distance_to(j) Measure the distance to j
equal: test if two iterators are equal
- Template Parameters:
derived_type –
value_t –
reference_t –
Public Types
-
using
self_type
= random_iterator_facade<derived_type, value_t, reference_t>¶
-
using
difference_type
= std::ptrdiff_t¶
-
using
reference
= reference_t¶
-
using
iterator_category
= std::random_access_iterator_tag¶
Public Functions
-
inline auto
operator[]
(const difference_type &n)¶
Friends
-
inline friend reference
operator*
(const derived_type &lhs)¶
-
inline friend const friend reference operator* (derived_type &lhs)
-
inline friend derived_type
operator++
(derived_type &lhs)¶
-
inline friend derived_type
operator++
(derived_type &lhs, int)¶
-
inline friend derived_type
operator--
(derived_type &lhs)¶
-
inline friend derived_type
operator--
(derived_type &lhs, int)¶
-
inline friend bool
operator==
(const derived_type &lhs, const derived_type &rhs)¶
-
inline friend bool
operator!=
(const derived_type &lhs, const derived_type &rhs)¶
-
inline friend bool
operator<=
(const derived_type &lhs, const derived_type &rhs)¶
-
inline friend bool
operator<
(const derived_type &lhs, const derived_type &rhs)¶
-
inline friend bool
operator>=
(const derived_type &lhs, const derived_type &rhs)¶
-
inline friend bool
operator>
(const derived_type &lhs, const derived_type &rhs)¶
-
inline friend auto
operator-
(const derived_type &lhs, const derived_type &rhs)¶
-
inline friend derived_type
operator+
(const derived_type &lhs, const difference_type &n)¶
-
inline friend derived_type
operator+
(const difference_type &n, const derived_type &rhs)¶
-
struct
hg::hg
::
region_adjacency_graph
¶ - #include <rag.hpp>
Result of the region adjacency graph (rag) construction algorithm.
-
template<typename
tree_t
, typenamenode_map_t
>
structhg::hg
::
remapped_tree
¶ - #include <common.hpp>
A simple structure to hold the result of a remapping operation on the nodes of a tree.
When algorithm transforms a tree into a new tree by removing or duplicating or reordering some of its nodes, it is useful to know the relation between the nodes of the new tree and the nodes of the original one.
For each node i of the new tree, node_map[i] gives the corresponding node in the original tree.
- Template Parameters:
tree_t –
node_map_t –
-
struct
hg::hg
::
scorer_partition_BCE
¶
-
struct
hg::hg
::
scorer_partition_DCovering
¶
-
struct
hg::hg
::
scorer_partition_DHamming
¶
-
template<typename
supervertex_labels_t
, typenametree_t
, typenamenode_map_t
>
structhg::hg
::
supervertex_hierarchy
¶ - #include <tree.hpp>
A simple structure to hold the result of supervertices_hierarchy algorithm.
- Template Parameters:
node_map_t –
tree_t –
node_map_t –
-
template<typename
transform_fun_type
, typenamebase_iterator_t
, typenamevalue_t
, typenamereference_t
= value_t>
structhg::hg
::
transform_forward_iterator
: public hg::forward_iterator_facade<transform_forward_iterator<transform_fun_type, base_iterator_t, value_t, value_t>, value_t, value_t>¶ - #include <iterators.hpp>
Facade to ease the declaration of new forward iterators as a transformation of an existing forward iterator (quick reimplementation of boost transform iterator).
- Template Parameters:
transform_fun_type –
base_iterator_t –
value_t –
reference_t –
Public Types
-
using
self_type
= transform_forward_iterator<transform_fun_type, base_iterator_t, value_t, reference_t>¶
Public Functions
-
inline void
increment
()¶
-
inline
transform_forward_iterator
(base_iterator_t &&base, transform_fun_type &&fun)¶
-
inline
transform_forward_iterator
(const base_iterator_t &base, const transform_fun_type &fun)¶
-
inline
transform_forward_iterator
()¶