Obtaining half-spaces from the convex hull in MATLAB

29 views (last 30 days)
Gayan Lankeshwara on 24 Feb 2022
Commented: Gayan Lankeshwara on 25 Feb 2022
I have determined the convex hull of certain points in 2-D plane using convhull function in MATLAB.
For example,
P = [0 0; 1 1; 1.5 0.5; 1.5 -0.5; 1.25 0.3; 1 0; 1.25 -0.3; 1 -1];
[k,av] = convhull(P);
will compute the 2-D convex hull of the points.
However, I am wondering whether there exists a direct approach (some kind of a function) to obtain the half-spaces determining the convex hull rather than manually calculating them?
Thank you.

Bruno Luong on 24 Feb 2022
Edited: Bruno Luong on 24 Feb 2022
% k = convhull(P), the half spaces form is {x such that A*x <= b} with
W=P(k,:);
V=W(1:end-1,:);
A=(W(2:end,:)-V)*[0 -1; 1 0]
b=dot(A,V,2)
Bruno Luong on 24 Feb 2022
Edited: Bruno Luong on 24 Feb 2022
P = [0 0; 1 1; 1.5 0.5; 1.5 -0.5; 1.25 0.3; 1 0; 1.25 -0.3; 1 -1];
k = convhull(P);
W=P(k,:);
V=W(1:end-1,:);
A=(W(2:end,:)-V)*[0 -1; 1 0];
b=dot(A,V,2);
% Check graphically, generate grid points
[xymi,xymax]=bounds(P,1);
x=linspace(xymi(1),xymax(1),33);
y=linspace(xymi(2),xymax(2),33);
[X,Y]=meshgrid(x,y);
XY=[X(:),Y(:)].';
% Check inside convhull using half-planes eqts
in=all(A*XY <= b,1);
figure
plot(W(:,1),W(:,2),'o-r');
hold on
plot(XY(1,in),XY(2,in),'.b');

John D'Errico on 24 Feb 2022
Edited: John D'Errico on 24 Feb 2022
What is a half space? Well, what do you mean by that? And what do you mean by "manually calculating" them? And is there some simple function you can call that will do exactly what you want? Probably not, but I don't really know what it is you want. No matter what, you would be using MATLAB, I presume. For example,
P = [0 0; 1 1; 1.5 0.5; 1.5 -0.5; 1.25 0.3; 1 0; 1.25 -0.3; 1 -1];
k = convhull(P)
k = 6×1
1 8 4 3 2 1
What does k mean here? k is an indirect index into the list of consecutive points on the convex hull. So the convex hull connexts points 1 and 8, then it moves to point 4, which implies a connection between points 4 and 8, etc.
P(k(1),:)
ans = 1×2
0 0
P(k(2),:)
ans = 1×2
1 -1
We can plot the points and the convex hull as...
plot(P(:,1),P(:,2),'bo',P(k,1),P(k,2),'r-')
If we subtract the pair of points, then we have a vector that is oriented along the direction of the vector along that edge. So the set of vectors that point along the edges for each of those line segments is simple to determine.
edgevecs = diff(P(k,:))
edgevecs = 5×2
1.0000 -1.0000 0.5000 0.5000 0 1.0000 -0.5000 0.5000 -1.0000 -1.0000
Better is to convert those vectors to each have unit length. So we might do
edgevecs = normalize(edgevecs,2,'norm')
edgevecs = 5×2
0.7071 -0.7071 0.7071 0.7071 0 1.0000 -0.7071 0.7071 -0.7071 -0.7071
That normalizes each row of the array of edgevecs to have unit length. But really, you may want normal vectors, assuming I think I know what you intend by the words half-spaces. Each line segment defines a line of infinite extent, and you want to know the half spaces that intersect to form the convex hull? So really, we want the normal vectors to those lines. This littel trick will give you the normal vectors:
edgenormals = edgevecs*[0 -1;1,0]
edgenormals = 5×2
-0.7071 -0.7071 0.7071 -0.7071 1.0000 0 0.7071 0.7071 -0.7071 0.7071
The matrix multiply I used there can be derived from a 2x2 rotation matrix, that describes a 90 deree rotation in the plane. It is a good trick to remember. Are those normal vectors inwards or outwards pointing normals? It is probably one or the other. :)
Given a point in the interior of the convex hull, we can compute a dot product to see which direction they point in.
C = mean(P(k(1:end-1),:),1)
C = 1×2
1 0
So the point C is essentially a point near the center of the convex hull. All that matters is C is INTERNAL to the convex hull, and by taking a convex linear combination of the nodes along that hull, we MUST have an internal point as a result. So which direction do those normals point along? Are they innies, or outies?
To determine that, we can just compute dot products with the edge normals, and the vectors pointing from C to a point on each of those edges. If that dot product is positive, then the vector was pointing outwards. If negative, then the vector was pointing towards the interior of the convex hull.
dotprod = sum((P(k(1:end-1),:) - C).*edgenormals,2)
dotprod = 5×1
0.7071 0.7071 0.5000 0.7071 0.7071
OK, so those vectors point outwards.
edgenormals = edgenormals.*sign(-dotprod)
edgenormals = 5×2
0.7071 0.7071 -0.7071 0.7071 -1.0000 0 -0.7071 -0.7071 0.7071 -0.7071
And now each of the vectors are inwards pointing normals. A point on each of those lines is simply
P(k(1:end-1),:)
ans = 5×2
0 0 1.0000 -1.0000 1.5000 -0.5000 1.5000 0.5000 1.0000 1.0000
And that completely defines the corresponding half-spaces, since a point along a line, and the normal vector to the line are sufficient.
Are you asking for a function that will now draw those domains in the plane, as colored regions? For that, you would probably want to write a little helper function that would take one normal vector and a point on the line, and then fill the region of interest in a figure. Then just call that function as many times as needed in a loop. That serves no purpose more than drawing a pretty picture of course. Is there a point to it?

Matt J on 24 Feb 2022
Edited: Matt J on 24 Feb 2022
You can use vert2lcon from
P = [0 0; 1 1; 1.5 0.5; 1.5 -0.5; 1.25 0.3; 1 0; 1.25 -0.3; 1 -1];
[A,b]=vert2lcon(P)
A = 5×2
1.0000 0.0000 0.7071 0.7071 0.7071 -0.7071 -0.7071 0.7071 -0.7071 -0.7071
b = 5×1
1.5000 1.4142 1.4142 0 0
The rows of A are the normals to the half spaces, pointing out of the hull. The equations for the half-spaces are A*x<=b.
Gayan Lankeshwara on 25 Feb 2022
Don't we need apply vert2lcon on the points forming the convex hull instead of P?

Categories

Find more on Bounding Regions in Help Center and File Exchange

R2019b

Community Treasure Hunt

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

Start Hunting!