Fast morphological reconstruction of large logical masks.

MATLAB's imreconstruct is slow for large 3D logical masks. bwreconstruct is a faster replacement.

MATLAB's imreconstruct() is optimised for 2D images (logical or numeric), but it suffers when working with large logical 3D (or higher dimension) matrices. bwreconstruct uses bwconncomp() to speed up this reconstruction considerably. It is around 10 times faster for one-off reconstructions, and faster still for multiple reconstructions (i.e., reconstructing a mask from a series of different markers).

IM = BWRECONSTRUCT(MARKER,MASK) performs morphological reconstruction
of the logical MARKER under the logical MASK. MARKER and MASK are two
two binary images with the same size, resulting in a binary image IM of
equivalent size.

IM = BWRECONSTRUCT(MARKER,MASKCC) performs similarly, but uses the
connected component structure MASKCC. MASKCC is the output from
BWCONNCOMP(MASK), and when it is already available, using MASKCC is
faster than providing a MASK binary image.

IM = IMRECONSTRUCT(MARKER,MASK,CONN) performs morphological
reconstruction with the specified connectivity. CONN may have the
following scalar values:

4 two-dimensional four-connected neighborhood
8 two-dimensional eight-connected neighborhood
6 three-dimensional six-connected neighborhood
18 three-dimensional 18-connected neighborhood
26 three-dimensional 26-connected neighborhood

Connectivity may be defined in a more general way (see IMRECONSTRUCT
for details).

Class support
-------------
MARKER and MASK must be logical arrays with the same dimension. If
MASKCC is used instead of MASK, it must relate to this same dimension.

Performance Note
----------------
This function performs faster than IMRECONSTRUCT for large binary
images of 3 or higher dimensions. For 2D or small images (sizes below
50-by-50-by-50, for example), IMRECONSTRUCT is still the faster option.

Example
---------
mask = false(100,100,100);
mask(:,[3 50 75],:) = true;
marker = false(size(mask));
marker(45:55,45:55,50:55) = true;
tic, im1 = imreconstruct(marker,mask); toc
tic, im2 = bwreconstruct(marker,mask); toc
isequal(im1, im2)
figure, patch(isosurface(im2,0.5),'FaceColor','g','EdgeColor','none'), hold on
patch(isosurface(marker,0.5),'FaceColor','b','EdgeColor','none')
patch(isosurface(mask & ~im2,0.5),'FaceColor','r','EdgeColor','none')
axis image, view(3), camlight

Note: I fully expect this to become a redundant function within a release or two of the IP toolbox!