Obtaining half-spaces from the convex hull in MATLAB

13 Ansichten (letzte 30 Tage)
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.

Akzeptierte Antwort

Bruno Luong
Bruno Luong am 24 Feb. 2022
Bearbeitet: Bruno Luong am 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)
  1 Kommentar
Bruno Luong
Bruno Luong am 24 Feb. 2022
Bearbeitet: Bruno Luong am 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');

Melden Sie sich an, um zu kommentieren.

Weitere Antworten (2)

John D'Errico
John D'Errico am 24 Feb. 2022
Bearbeitet: John D'Errico am 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
Matt J am 24 Feb. 2022
Bearbeitet: Matt J am 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.
  1 Kommentar
Gayan Lankeshwara
Gayan Lankeshwara am 25 Feb. 2022
Don't we need apply vert2lcon on the points forming the convex hull instead of P?

Melden Sie sich an, um zu kommentieren.

Kategorien

Mehr zu Bounding Regions finden Sie in Help Center und File Exchange

Produkte


Version

R2019b

Community Treasure Hunt

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

Start Hunting!

Translated by