diff --git a/src/algorithms.hpp b/src/algorithms.hpp
index 797f471bc51afac58f650d728ea57cfec445e4ab..39df831c26a8c8a0321128a19396ea6a5d092b30 100644
--- a/src/algorithms.hpp
+++ b/src/algorithms.hpp
@@ -112,29 +112,86 @@ namespace algorithms{
}
template<typename C>
- bool is_contiguously_numbered(const C &parent_list)
+ bool has_contiguous_segments(const C &parent_index)
{
static_assert(
std::is_integral<typename C::value_type>::value,
"integral type required"
);
- std::vector<bool> is_leaf(parent_list.size(), false);
+ if (!is_minimal_degree(parent_index)) {
+ return false;
+ }
+
+ std::vector<bool> is_leaf(parent_index.size(), false);
- auto ret = true;
- for (std::size_t i = 1; i < parent_list.size(); ++i) {
- if (is_leaf[parent_list[i]]) {
- ret = false;
- break;
+ for (std::size_t i = 1; i < parent_index.size(); ++i) {
+ auto p = parent_index[i];
+ if (is_leaf[p]) {
+ return false;
}
- if (parent_list[i] != i-1) {
+ if (p != i-1) {
// we have a branch and i-1 is a leaf node
is_leaf[i-1] = true;
}
}
- return ret;
+ return true;
+ }
+
+ template<typename C>
+ std::vector<typename C::value_type> child_count(const C &parent_index)
+ {
+ static_assert(
+ std::is_integral<typename C::value_type>::value,
+ "integral type required"
+ );
+
+ std::vector<typename C::value_type> count(parent_index.size(), 0);
+ for (std::size_t i = 1; i < parent_index.size(); ++i) {
+ ++count[parent_index[i]];
+ }
+
+ return count;
+ }
+
+ template<typename C, bool CheckStrict = true>
+ std::vector<typename C::value_type> branches(const C &parent_index)
+ {
+ static_assert(
+ std::is_integral<typename C::value_type>::value,
+ "integral type required"
+ );
+
+ if (CheckStrict && !has_contiguous_segments(parent_index)) {
+ throw std::invalid_argument(
+ "parent_index has not contiguous branch numbering"
+ );
+ }
+
+ auto num_child = child_count(parent_index);
+ std::vector<typename C::value_type> branch_index(
+ parent_index.size(), 0
+ );
+
+ std::size_t num_branches = (num_child[0] == 1) ? 1 : 0;
+ for (std::size_t i = 1; i < parent_index.size(); ++i) {
+ auto p = parent_index[i];
+ if (num_child[p] > 1) {
+ ++num_branches;
+ }
+
+ branch_index[i] = num_branches;
+ }
+
+ return branch_index;
+ }
+
+ template<typename C>
+ std::vector<typename C::value_type> branches_fast(const C &parent_index)
+ {
+ return branches<C, false>(parent_index);
}
} // namespace algorithms
diff --git a/tests/test_algorithms.cpp b/tests/test_algorithms.cpp
index 91e43c8b86717ce536110caa7394b9b8b6d74d0d..05da8c0c41f8aeb39c29f35f65c814542cd1355b 100644
--- a/tests/test_algorithms.cpp
+++ b/tests/test_algorithms.cpp
@@ -170,7 +170,7 @@ TEST(algorithms, is_positive)
);
}
-TEST(algorithms, is_contiguously_numbered)
+TEST(algorithms, has_contiguous_segments)
{
//
// 0
@@ -178,13 +178,13 @@ TEST(algorithms, is_contiguously_numbered)
// 1
// |
// 2
- // /|\
+ // /|\.
// 3 7 4
- // / \
+ // / \.
// 5 6
//
EXPECT_FALSE(
- nest::mc::algorithms::is_contiguously_numbered(
+ nest::mc::algorithms::has_contiguous_segments(
std::vector<int>{0, 0, 1, 2, 2, 3, 4, 2}
)
);
@@ -201,7 +201,7 @@ TEST(algorithms, is_contiguously_numbered)
// 4 7
//
EXPECT_FALSE(
- nest::mc::algorithms::is_contiguously_numbered(
+ nest::mc::algorithms::has_contiguous_segments(
std::vector<int>{0, 0, 1, 2, 3, 2, 2, 5}
)
);
@@ -218,7 +218,7 @@ TEST(algorithms, is_contiguously_numbered)
// 4 6
//
EXPECT_TRUE(
- nest::mc::algorithms::is_contiguously_numbered(
+ nest::mc::algorithms::has_contiguous_segments(
std::vector<int>{0, 0, 1, 2, 3, 2, 5, 2}
)
);
@@ -235,23 +235,130 @@ TEST(algorithms, is_contiguously_numbered)
// 4 6
//
EXPECT_TRUE(
- nest::mc::algorithms::is_contiguously_numbered(
+ nest::mc::algorithms::has_contiguous_segments(
std::vector<int>{0, 0, 1, 2, 3, 2, 5, 1}
)
);
// Soma-only list
EXPECT_TRUE(
- nest::mc::algorithms::is_contiguously_numbered(
+ nest::mc::algorithms::has_contiguous_segments(
std::vector<int>{0}
)
);
// Empty list
EXPECT_TRUE(
- nest::mc::algorithms::is_contiguously_numbered(
+ nest::mc::algorithms::has_contiguous_segments(
std::vector<int>{}
)
);
+}
+
+TEST(algorithms, child_count)
+{
+ {
+ //
+ // 0
+ // /|\
+ // 1 4 6
+ // / | \
+ // 2 5 7
+ // / \
+ // 3 8
+ // / \
+ // 9 11
+ // / \
+ // 10 12
+ // \
+ // 13
+ //
+ std::vector<int> parent_index =
+ { 0, 0, 1, 2, 0, 4, 0, 6, 7, 8, 9, 8, 11, 12 };
+ std::vector<int> expected_child_count =
+ { 3, 1, 1, 0, 1, 0, 1, 1, 2, 1, 0, 1, 1, 0 };
+
+ // auto count = nest::mc::algorithms::child_count(parent_index);
+ EXPECT_EQ(expected_child_count,
+ nest::mc::algorithms::child_count(parent_index));
+ }
+
+}
+
+TEST(algorithms, branches)
+{
+ using namespace nest::mc;
+
+ {
+ //
+ // 0
+ // /|\
+ // 1 4 6
+ // / | \
+ // 2 5 7
+ // / \
+ // 3 8
+ // / \
+ // 9 11
+ // / \
+ // 10 12
+ // \
+ // 13
+ //
+ std::vector<int> parent_index =
+ { 0, 0, 1, 2, 0, 4, 0, 6, 7, 8, 9, 8, 11, 12 };
+ std::vector<int> expected_branches =
+ { 0, 1, 1, 1, 2, 2, 3, 3, 3, 4, 4, 5, 5, 5 };
+ auto actual_branches = algorithms::branches_fast(parent_index);
+ EXPECT_EQ(expected_branches, actual_branches);
+ }
+
+ {
+ //
+ // 0
+ // |
+ // 1
+ // |
+ // 2
+ // |
+ // 3
+ //
+ std::vector<int> parent_index =
+ { 0, 0, 1, 2 };
+ std::vector<int> expected_branches =
+ { 0, 1, 1, 1 };
+
+ auto actual_branches = algorithms::branches_fast(parent_index);
+ EXPECT_EQ(expected_branches, actual_branches);
+ }
+
+ {
+ //
+ // 0
+ // |
+ // 1
+ // |
+ // 2
+ // / \
+ // 3 4
+ // \
+ // 5
+ //
+ std::vector<int> parent_index =
+ { 0, 0, 1, 2, 2, 4 };
+ std::vector<int> expected_branches =
+ { 0, 1, 1, 2, 3, 3 };
+
+ auto actual_branches = algorithms::branches_fast(parent_index);
+ EXPECT_EQ(expected_branches, actual_branches);
+ }
+
+ {
+ std::vector<int> parent_index = { 0 };
+ std::vector<int> expected_branches = { 0 };
+
+ auto actual_branches = algorithms::branches_fast(parent_index);
+ EXPECT_EQ(expected_branches, actual_branches);
+ }
}