Shouldn't bwconvhull() be idempotent?

4 views (last 30 days)
Matt J
Matt J on 20 Jun 2021
Edited: Matt J on 1 Jul 2021
Shouldn't the following two calculations give the same results?
load BWimage
A=bwconvhull(BW);
B=bwconvhull(bwconvhull(BW));
They don't.
isequal(A,B)
ans = logical
0

Accepted Answer

Steve Eddins
Steve Eddins on 30 Jun 2021
Hi Matt,
At a general level - and speaking from too much personal experience - expectations about geometrical operations such as the convex hull are often not met when we're dealing with a discrete grid instead of a continuous domain. Since bwconvhull is inherently working on inputs and outputs on a discrete grid, it is really only an approximation of a true convex hull, and I wouldn't expect it to be exactly idempotent.
With approximations, different approximation techniques usually have their own pros and cons, and one can consider whether one approximation technique is better than another, and under what conditions.
Your function treats pixels as points. It computes the convex hull of all the pixel centers and then uses inpolygon to "rasterize" the resulting shape. A disadvantage of that approach is that it will make skinny, single-pixel lines disappear competely because the convex hull of the pixel centers is a degenerate polygon with no area. On the following input, your function produces no foreground pixels at all in the result.
W = zeros(30,30);
W(15,10:20) = 1;
I don't think we could ship this method without somehow addressing that issue.
The function bwconvhull treats pixels as squares having unit area. It computes the convex hull of the set of all the pixel corners and then uses roipoly to rasterize the result. This method behaves better for skinny shapes, but it behaves less well in terms of idempotence, as you have clearly demonstrated.
Oversimplifying a bit, it comes down to tie-breaking rules related to polygons that go exactly through a pixel center. Is such a pixel inside or outside the polygon? Each choice results in its own set of nonideal behaviors.
It is worth asking, though, whether the method in bwconvhull can be improved with respect to idempotence while avoiding the degeneracy I mentioned above. We can take a look at this.
  5 Comments
Steve Eddins
Steve Eddins on 1 Jul 2021
Thanks, Matt. I will take a further look at your suggestions to see if there is something useful we can do.
Matt J
Matt J on 1 Jul 2021
Thank you, too. It was a good discussion.

Sign in to comment.

More Answers (2)

Jonas
Jonas on 20 Jun 2021
Edited: Jonas on 20 Jun 2021
ideally it should. already the shape you provide in the original is convex, but the result of the first call of bwconvhull and the original differ at 130 positions. the result of the first bwconvhull result and second show differences at 114 positions.
i think the general issue is that shapes/lines can not be represented good in pixel format except for horizontal, vertical and diagonal lines. in the other cases the algorithm just rounds some pixel to be inside the hull or not
using bwconvhull on easier shapes like a block or something like
[0 0 0 0;
0 1 1 0;
0 1 1 1;
0 0 1 0];
is idempotent
  2 Comments
Matt J
Matt J on 20 Jun 2021
Edited: Matt J on 20 Jun 2021
i think the general issue is that shapes/lines can not be represented good in pixel format except for horizontal, vertical and diagonal lines.
Sure, but the test below shows that whatever method bwconvhull is using, it is strangely and unnecessarily unstable. WIth my own implementation, I am able to iterate the convex hull operation many times without distorting the original shape nearly so much as bwconvhull does.
load tst
B=A;
tic;
for i=1:100; B=bwconvhull(B); end
toc;
Elapsed time is 0.475528 seconds.
C=A;
tic;
for i=1:100; C=bwconvhullMatt(C); end
toc
Elapsed time is 2.404097 seconds.
tiledlayout(1,3,'Padding','none');
nexttile, imshow(A); title 'Original'
nexttile, imshow(B); title 'BWCONVHULL'
nexttile, imshow(C); title 'Matt''s BWCONVHULL'
function Ch=bwconvhullMatt(BW)
persistent I J
if isempty(I)
sz=size(BW);
[I,J]=ndgrid(1:sz(1),1:sz(2));
end
V=cell2mat(bwboundaries(BW));
warning off
pgon=convhull(polyshape(V));
warning on
Ch=reshape( inpolygon(I(:),J(:),pgon.Vertices(:,1), pgon.Vertices(:,2)), size(I) );
end
Jonas
Jonas on 20 Jun 2021
another example where built in matlab functions work worse than user developed functions.

Sign in to comment.


Matt J
Matt J on 1 Jul 2021
Edited: Matt J on 1 Jul 2021
Below is a 3rd version of the function I proposed in my various conversations with Steve and Jonas. Assuming one is happy with treating the convex hull as the hull of the pixel centers, this version will compute it faster than the previous versions, and without the problems handling skinny 1 pixel lines. It does, however, have the issue with adjacent regions mentioned by Steve, though an imclose() operation could be used to fix that.
load tst
B=A;
tic;
for i=1:100; B=bwconvhull(B); end
toc;
Elapsed time is 0.690813 seconds.
C=A;
tic;
for i=1:100; C=bwconvhullCenters(C); end
toc
Elapsed time is 0.809547 seconds.
tiledlayout(1,3,'Padding','none');
nexttile, imshow(A); title 'Original'
nexttile, imshow(B); title 'bwconvhull'
nexttile, imshow(C); title 'bwconvhullCenters'
function Ch=bwconvhullCenters(BW)
d=0.0001;
dr=[-d,-d,+d,+d];
dc=[-d,+d,-d,+d];
V=bwboundaries(BW);
[r,c]=deal(V{1}(:,1), V{1}(:,2));
r=r+dr;
c=c+dc;
k=convhull(r,c,'Simplify',true);
r=r(k); c=c(k);
Ch=roipoly(BW,c,r);
end

Products


Release

R2020b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!