VoxelComplex: CriticalKernels framework

[phcerdan]

PR Description

Implementation of Critical Kernel framework based on M.Couprie and G.Bertrand articles:

• Asymmetric parallel 3D thinning scheme and algorithms based on isthmuses.
• Powerful Parallel and Symmetric 3D Thinning Schemes Based on Critical Kernels
• A 3D Sequential Thinning Scheme Based on Critical Kernels
• Add general description of framework for docs.
• Add LookUpTable interface from Add LookUpTableFunctions #1155 (speeds an order of magnitude)

CLI script for testing in any image (future PR in DGtalTools)

Implementation Details

• VoxelComplex inherits from CubicalComplex, but no override or hide any method.

• It has an Object member storing the voxels(spels) of an input image. This allows the use of the isSimple method from Object, and also the ObjectBoostGraphInterface.

• CriticalCliques:
  A clique is a set of Voxels.
  K_0, K_1, K_2, and K_3 return a pair<bool, Clique> where Clique=CubicalComplex, and the bool is true if the Clique is critical (as described in papers).

  K_3 applies to Voxels, and it uses the isSimple method of Object.

  K_2 applies to surfels, K_1 to linels, and K_0 to pointels.

In VoxelComplexFunctions a couple of thinning methods are implemented in a functional style, taking a VoxelComplex, a Select function, and a Skel function.

Select function: Select a voxel from a clique (set of voxels)

• Random.
• First: just the first voxel of the set.
• MaxValue: Wrap to use with a lambda, to use a DGtal::DistanceTransformation to select the voxel with greatest distance value.

Skel function: Predicate to check if input voxel is part of the skeleton that we are interested to preserve.

• Ultimate: Null, preserve only to keep same topology. The result of the thinning in a connected object is just one voxel.
• End: preserve if voxel has only one neighbor.
• OneIsthmus: preserve if voxel is oneIsthmus, as described in articles.
• TwoIsthmus: idem, this conserves planar structures.
• IsIsthmus : (OneIsthmus || TwoIsthmus)
• withTable: Wrap to use with a lambda, that accept any table from LookUpTable ( Add LookUpTableFunctions #1155 ). Test provides example of how to use with isIsthmus table.

WIP Todo list:

• In persistence algorithm.
  • [NA] provide reliable test for persistence. Hard to test results on this. The algorithm is simple enough, and the results make sense.
• Optimize:
  • [NA] Ensure that openMP code in criticalCliques works in parallel in clang (no compile errors so far). Left, with LUT criticalCliques speed is not that important.
• Solve all TODO. :vimgrep TODO .ih src/*
  • General Method to get spels from Cells. SOLVED: Kneighborhood. remove comment.
• Remove lexicographical order from Select functions.

  Checklist

• Unit-test of your feature with Catch.
• Doxygen documentation of the code completed (classes, methods, types, members...)
  • Revise criticalKernel documentation, specially interface to K_2.
• Documentation module page added or updated.
• New entry in the ChangeLog.md added.
• No warning raised in Debug cmake mode (otherwise, Travis C.I. will fail).
• All continuous integration tests pass (Travis & appveyor)

[JacquesOlivierLachaud]

Excellent idea to implement the works of Bertrand and Couprie. I do not have much time at the moment to look at your PR, but I will definitely do it later. By the way, Object::isSimple should also be tabulated (in fact, there is an example that builds the LUT).

[phcerdan]

@JacquesOlivierLachaud What do you think about distributing in the source code those pre-computed LUT plain files?
And for the user to load those tables... Should we read the file with a ifstream or with boost::serialization::bitset on the LUT's?

Just mentioning boost::serialization for the binary_file option, which is faster to read than a plain file. But probably riskier in terms of portability. Another option is to generate the binary file at build time from the distributed plain file...

Edited: I will test the options when the LUT generation finishes in my computer, maybe it is fast enough to read the table in plain text...

[phcerdan]

Look Up Tables ( #1155 ) of simplicity and isIsthmus speeds up the computation a lot.
I have merged #1155 in this PR.

[phcerdan]

Appveyor test looks as a failure because some timeout. I will try to solve all the checklist soon.

[phcerdan]

TODO: Change tests to use LUT.
AppVeyor fails because timeout. Tests are not using look up tables yet, and they are taking a while.

[phcerdan]

Hey, I will come back to this pretty soon to polish it. I am not sure if it will ok for dealing with large images, but I will give it a try.

[dcoeurjo]

👌

[phcerdan]

@rolanddenis I started this before some changes lead by you with khalimsky spaces. I adapted my code to compile against the changes, but not sure if I am taking full advantage of those improvements. Any hint is appreciated! 😅

The code I think should be ready for a preliminary review. The documentation of this PR is pretty much enough for its usage, please see the associated branch in DGtalTools for a screenshot (it needs this PR to work).

The script in DGtalTools criticalKernelsThinning3D.

The code is quite flexible, accepting a skel function for keeping the voxels that the user want to keep in the skeletonization process.

THe idea of isthmus is based on Couprie and Bertrand work, and it is 3D only, but it seems that has been ported to 2D as well here: 2D parallel thinning algorithms based on isthmus-preservation, but I have not pursued any work in this direction.

Compared with the homotopicThinning, this code is able to get more details (the finger are captured). Also it is able to use the persistence parameter as a pretty consistent trimmer of short/noise branches.

Any feedback is welcome, I will be glad to push this in the next few days if our time permits. Cheers!

[dcoeurjo]

hi @phcerdan can this PR be reviewed?

[phcerdan]

Yes, it can be!

[rolanddenis]

Sorry @phcerdan, I just notice your question about the changes in Khalimsky space ...

To be short, the main difference between KhalimskySpaceNDand KhalimskyPreSpaceND is that the former has bounds (i.e. a domain) and the latter is unbounded. The KhalimskyPreSpaceND was also introduced to allow direct manipulation of cell coordinates without risking of having "invalid" cells. As a consequence, KhalimskySpaceND has been modified to prohibit such usages and to guarantee that cells are always valid (at least in debug mode).

In the class VoxelComplexyou introduced here, I think there is yet no problem since it inherits from CubicalComplexwhich needs a Khalimsky space that respects the CCellularGridSpaceND concept (i.e. only KhalimskySpaceND).

Maybe CubicalComplex could be modified to accept model of CPreCellularGridSpaceND instead (what do you think @JacquesOlivierLachaud ?) but DGtal::Object also requires a domain : VoxelComplex::populateObjectFromCells actually use CCellularGridSpaceNDspecific methods (lowerBound and upperBound) to instantiate a TObject.

Finally, while fast checking the code, I didn't see any wrong usage of KhalimskySpaceND 😉

[rolanddenis]

What I do sometimes in this case is to propose providing a whole random generator instead of the seed only, something like that :

template < typename TComplex >
std::pair<typename TComplex::Cell, typename TComplex::Data>
DGtal::functions::selectRandom(
      const typename TComplex::Clique & clique)
{
    static std::mt19937 gen( std::random_device{}() );
    return selectRandom(clique, gen);
}

template < typename TComplex, typename TRandomGenerator >
std::pair<typename TComplex::Cell, typename TComplex::Data>
DGtal::functions::selectRandom(
      const typename TComplex::Clique & clique,
      TRandomGenerator & gen
)
{
    auto size = clique.nbCells(3);
    std::uniform_int_distribution<> dis(0, size - 1);
    auto it = clique.begin(3);
    std::advance(it, dis(gen));
    return *it;
}

[phcerdan]

Good idea, changing it.

[rolanddenis]

Oops, I forget to specify the template in my suggestion :

return selectRandom<TComplex>(clique, gen);

[dcoeurjo]

Could you please keep me updated with this PR (Release pending) ? Thanks

[phcerdan]

@dcoeurjo I think is ready to be merged.
The docs could be improved with some images from this DGtalTools PR once is merged.

[dcoeurjo]

👌 thanks

[dcoeurjo]

Would you please add further comment on the thinning algorithms ?

[dcoeurjo]

hi @phcerdan
I'm about to release the lib. Could you please update the doc ? thanks

[phcerdan]

Hi @dcoeurjo, sorry for the delay on this. I will do it first thing in the morning, one more day please!

[phcerdan]

@dcoeurjo docs updated! Thanks for the patience and the review. I think it is ready.

[phcerdan]

Errors in travis seem unrelated. Time out and doxygen not found.

[dcoeurjo]

thanks a lot!

[dcoeurjo]

allgood;)
thanks a lot

[phcerdan]

🎉