Question about gmmTreeRegEstep function in HGMM Reg about probreg

[SpriteLin-ZJU]

I think the following changes should be correct? Executing accumulate function before break.

accumulate(moments_private[search_id], gamma[search_id - j0], points.row(i)); if (complexity(std::get<2>(nodes[search_id])) <= lambda_c) break;

[neka-nat]

Hi @SpriteLin-ZJU ,

Thanks you for the report!
You're right, the implementation is wrong.
I think the following implementation is correct.

    for (Integer i = 0; i < points.rows(); ++i) {
        Integer search_id = -1;
        Vector gamma = Vector::Zero(N_NODE);
        Integer j0 = 0;
        for (Integer l = 0; l < max_tree_level; ++l) {
            j0 = child(search_id);
            for (Integer j = j0; j < j0 + N_NODE; ++j) {
                gamma[j - j0] = std::get<0>(nodes[j]) *
                                gaussianPdf(points.row(i), std::get<1>(nodes[j]), std::get<2>(nodes[j]));
            }
            const Float den = gamma.sum();
            if (den > eps) {
                gamma /= den;
            } else {
                gamma.fill(0.0);
            }
            gamma.maxCoeff(&search_id);
            search_id += j0;
            if (complexity(std::get<2>(nodes[search_id])) <= lambda_c) break;
        }
        accumulate(moments[search_id], gamma[search_id - j0], points.row(i));
    }

I also confirmed that it works.

[SpriteLin-ZJU]

@neka-nat Thanks for reply.
I think it workd now.

But I am still puzzled by the registration E step in the paper. Why does the algorithm only calculate the moments of the leaf node, but not the moments of the parent node?

[neka-nat]

My guess is that updating just the leaf node or the whole thing would make little difference to the results.
Since the parent's node information is the sum of the children node information, it would be no different than using only the child's information.

[SpriteLin-ZJU]

@neka-nat
I think you are right.
I still have a few questions about the code implementation:

1. In the buildGmmTree(const MatrixX3& points):

   auto idxs = (n_total * Vector::Random(n_total)).array().abs().cast();
   Maybe you mean
   auto idxs = (points.rows() * Vector::Random(n_total)).array().abs().cast();

2. For a large point cloud (200k), calculating the sig2 is too time-consuming.

3. Also in buildGmmTree:
   The variable q will increase as the number of points in the point cloud increases, and may eventually increase to several orders of magnitude more than lambda_s. As a result, the convergence speed is slow or even unable to converge when building trees. It may be more effective to choose a variable that can dynamically change with the number of points instead of a fixed value.

[neka-nat]

Thank you for reporting.
I've fixed each bug and pushed to master.

1. n_total -> points.rows()
2. I removed the points.rows() loop and initialize the covariance matrix in the n_total loop.
3. I added lamda_s parameter.

[SpriteLin-ZJU]

Please check if the second modification is correct, it gives the wrong result in my test data.

Here is my test data.
waymo.pcd.zip

[SpriteLin-ZJU]

I used the following code to initialize the covariance matrix. But it seems that the randomness of the results has become greater, and sometimes the registration fails. I don't know what caused this.

#pragma omp parallel for num_threads(THREAD_NUM)
for (Integer j = 0; j < n_total; ++j)
{
std::get<0>(nodes[j]) = 1.0 / N_NODE;
std::get<1>(nodes[j]) = points.row(idxs[j]);
MatrixX3 tmp = (points.rowwise() - points.row(idxs[j])).matrix();
std::get<2>(nodes[j]) = tmp.transpose() * tmp / points.rows();
}

[neka-nat]

The initialization has been improved to take the tree structure into account.
The results still depend on the initial values, but can be stabilized by adjusting parameters, tree_level and maxiter.
Please try the latest master.

[SpriteLin-ZJU]

Hi, @neka-nat ,
I haven't tested the effect of different initialization methods on the registration effect, but I think the following modifications should be made in the nodes initialization function.

1.for (Integer j = 0; j < max_tree_level * N_NODE; ++j)
------->
for (Integer j = 0; j < std::pow(N_NODE, max_tree_level); ++j)

2.for (Integer j = 0; j < N_NODE * (l + 1); ++j)
------->
for (Integer j = 0; j < std::pow(N_NODE, l+1); ++j)

[neka-nat]

Thanks for the correction.
I've pushed the fix to master.
The results have been more stable for the stanford bunny.