Region splitting

[annimesh2809]

Added functions for creating region trees (using RegionTrees.jl) and segmenting images using region splitting algorithm. The time complexity of the algorithm is O(n*log(n)). It works for N-dimensional arbitrary indexed images.

Working example:

using Images, ImageSegmentation, TestImages
img = testimage("camera");
f{T}(a::T) = Images.accum(T)(a);
homogeneous(i::AbstractArray)::Bool = (abs(f(minimum(i))-f(maximum(i))) < 0.2);
seg = region_splitting(img, homogeneous);

Original:

RegionTree (threshold = 0.2):

Segmented Image:

[codecov[bot]]

Codecov Report

Merging #10 into master will decrease coverage by 0.4%. The diff coverage is 95.55%.

@@            Coverage Diff             @@
##           master      #10      +/-   ##
==========================================
- Coverage    99.5%   99.09%   -0.41%     
==========================================
  Files           6        7       +1     
  Lines         400      444      +44     
==========================================
+ Hits          398      440      +42     
- Misses          2        4       +2

Impacted Files	Coverage Δ
src/ImageSegmentation.jl	100% <ø> (ø)	:arrow_up:
src/core.jl	100% <100%> (ø)	:arrow_up:
src/region_merging.jl	95.45% <95.45%> (ø)

---

Continue to review full report at Codecov.

Legend - Click here to learn more Δ = absolute <relative> (impact), ø = not affected, ? = missing data Powered by Codecov. Last update 02ebd26...cd415b6. Read the comment docs.

[timholy]

:+1: Does the use of the MVector make those splats inferrable? Assuming yes, it's a nice trick.

I didn't see any particular issues, other than the obvious thought that there's a bit of code duplication. But it's not a very long stretch, so I'm fine with it.

[annimesh2809]

I can't understand why there are so many red bars in ProfileView although @code_warntype doesn't give any warnings.

Code used for profiling:

using Images, ImageSegmentation, ProfileView, TestImages
img = testimage("lena_gray_512")
f{T}(a::T) = Images.accum(T)(a);
homogeneous(i::AbstractArray)::Bool = (abs(f(minimum(i))-f(maximum(i))) < 0.2);
t = region_tree(img, homogeneous)
Profile.clear()
@profile region_tree(img, homogeneous)
ProfileView.view()

Output:

I tried to debug it but could not understand the origin of the problem. There is also a decrease in performance as compared to the "only 2-d implementation". Please help.

[timholy]

So, region_tree itself (and all the code it inlines) is inferrable, but region_tree! is not:

julia> rtree = Cell(SVector(first(CartesianRange(indices(img))).I), SVector(length.(indices(img))), (0.,0))
Cell: RegionTrees.HyperRectangle{2,Float64}([1.0, 1.0], [512.0, 512.0])

julia> @code_warntype ImageSegmentation.region_tree!(rtree, img, homogeneous)
...

So, the problem seems to be that the MVector trick doesn't actually make the call inferrable (even though in principle it could). You could add something like

rvt = tuple(rv...)::NTuple{N,UnitRange{Int}}
imgv = view(img, rvt...)

but you could also create bv and rt as tuples directly using "lispy tuple programming," which makes it quite easy for the compiler to reason about the types of objects. I'll illustrate how you might do this for creating bv, and leave rt as an exercise to make sure you have a chance to understand this:

julia> start_ind = CartesianIndex(1,2,3,4,5)  # this will encode the length of the tuple we want
CartesianIndex{5}((1, 2, 3, 4, 5))

julia> makebt(i) = ()     # this will be called when we have processed the final index (when the splat below is empty)
makebt (generic function with 1 method)

julia> @inline makebt(i, ind, rest...) = (i & 1 != 0, makebt(i>>1, rest...)...)  # forced inlining means this effectively runs at *compile time*, not runtime
makebt (generic function with 2 methods)

julia> makebt(27, start_ind.I...)
(true, true, false, true, true)

julia>  bits(27)
"0000000000000000000000000000000000000000000000000000000000011011"

Also check out the @code_llvm for this call, you'll see it's impressively compact and allocates only for the return value (which gets elided in the actual compiled code).

If you have any questions, don't hesitate to ask: lispy programming is a different way of thinking, but it's extremely powerful (essentially all of Julia's internal array-handling code is written this way) and when used correctly delivers awesome performance: it gets rid of all looping/branching and acts like you've hand-coded each dimensionality for optimal performance.

@tejus-gupta, it would be worth making sure you understand this, too.

[timholy]

Also worth pointing out if you don't realize this: you can "create" an infinite number of tuples (of known type/length) without ever allocating memory, because the compiler elides them. (The fact that you don't have to allocate memory on each iteration is the only reason CartesianIndex can work.) This is very handy for many indexing operations.

[annimesh2809]

:tada: That's a really nice way of programming (and powerful too). It fixed the major source of instability. The performance has increased 10 folds.

There are still some minor instability issues related to mean and mapfoldl. Let me try to debug them, if I fail then I will ask for help. :smile:

[timholy]

Glad that worked!

Not relevant for you now, but just FYI in 0.7/1.0 you'll want to omit the @inline because the inliner is now smarter (and leaving it out can actually solve some problems with excessive compile time at high dimensionality, e.g., https://github.com/JuliaLang/julia/pull/22989#discussion_r129979517).

[codecov-io]

Codecov Report

Merging #10 into master will decrease coverage by 0.4%. The diff coverage is 95.55%.

@@            Coverage Diff             @@
##           master      #10      +/-   ##
==========================================
- Coverage    99.5%   99.09%   -0.41%     
==========================================
  Files           6        7       +1     
  Lines         400      444      +44     
==========================================
+ Hits          398      440      +42     
- Misses          2        4       +2

Impacted Files	Coverage Δ
src/ImageSegmentation.jl	100% <ø> (ø)	:arrow_up:
src/core.jl	100% <100%> (ø)	:arrow_up:
src/region_merging.jl	95.45% <95.45%> (ø)

---

Continue to review full report at Codecov.

Legend - Click here to learn more Δ = absolute <relative> (impact), ø = not affected, ? = missing data Powered by Codecov. Last update 02ebd26...15cff28. Read the comment docs.

[timholy]

Woot! Thanks for this!