I don't have access to the code, so I can only surmise (since I've written alpha shape codes before in MATLAB, but not the one in question) ...
From the help, we see:
SHP = alphaShape(X,Y,A) creates a 2D alpha shape of the points (X,Y),
with alpha radius A. X and Y are column vectors that specify the
coordinates. A is a positive scalar that defines the alpha radius. If
A is not specified, the default is the smallest alpha that produces an
alpha shape enclosing all of the points.
The basic idea behind an alpha shape is that we start with a delaunay triangulation. We can think of that as an alpha shape with an alpha ball of infinite radius
xy = randn(100,2);
>> S = alphaShape(xy,inf)
S =
alphaShape with properties:
Points: [100×2 double]
Alpha: Inf
HoleThreshold: 0
RegionThreshold: 0
The alpha ball erodes into the triangulation. If a ball of radius alpha can penetrate a triangle and touch an interior point, then the triangle is effectively deleted from the list. (It has been many years since I wrote alpha shape codes of my own, but this is the essential idea.)
We can see however, that if we reduce alpha from infinity, that gradually we will eliminate those triangles around the edges of the domain. That is the goal of an alpha shape. But how far do we go? The default, as we are told, is the radius that does not cause any points to be deleted from the shape.
S = alphaShape(xy)
S =
alphaShape with properties:
Points: [100×2 double]
Alpha: 0.59772
HoleThreshold: 0
RegionThreshold: 0
plot(S)
hold on
plot(xy(:,1),xy(:,2),'ro')
So the default value was chosen to be 0.59772. Smaller than that? The alpha ball can erode some points effectively out of existence, by deleting all triangles tht contined those points s vertices.
S2 = alphaShape(xy,0.4)
S2 =
alphaShape with properties:
Points: [100×2 double]
Alpha: 0.4
HoleThreshold: 0
RegionThreshold: 0
plot(S2)
hold on
plot(xy(:,1),xy(:,2),'ro')
However, again, I don't have the code. And it looks like the choice of a default for alpha is a moderately crude one. Since if I choose a slightly smaller value for alpha than the default, the shape still contains all the original points as vertices.
S3 = alphaShape(xy,0.58)
S3 =
alphaShape with properties:
Points: [100×2 double]
Alpha: 0.58
HoleThreshold: 0
RegionThreshold: 0
But if you plot it, visually, this seems to be identical to that from the default value. That means the phrase we read:
the default is the smallest alpha that produces an
alpha shape enclosing all of the points
was only an approximation.
So clearly the scheme used to compute a default alpha is a moderately crude one, ad hoc, but not an optimization of any sort. I can imagine some schemes I might choose. What I need to recognize is the author of the code wanted to choose a scheme that was efficient. That suggests a sort of greedy algorithm. Once an alpha is identified that does reasonably well by some criteria, the search stops.
The question might be, what would I have done? I might start with the list of edges.
tri = DT.ConnectivityList;
edgelist = [tri(:,[1 2]);tri(:,[1 3]);tri(:,[2 3])];
edgelist = unique(sort(edgelist,2),'rows');
edgelens = sqrt(sum(xy(edgelist(:,1),:) - xy(edgelist(:,2),:)).^2,2);
Now, get the median edgelength.
median(edgelens)
ans =
0.36336
So this gives me some measure of what a number that is surely going to be too small for alpha. I might try creating an alpha shape with that as alpha. If no points get deleted, then I would be done. Otherwise, I might just increase that value by some fraction, testing to look for a value where no points get deleted. But clearly the search that is done is an ad hoc, approximate one. And sadly, I cannot know better than that. (Since I am not an employee.)
Best Answer