geometry icon indicating copy to clipboard operation
geometry copied to clipboard

Point to hull distance

Open Cedvega opened this issue 4 years ago • 8 comments

Would it be possible to add the point to hull distance ?

Cedvega avatar Aug 25 '21 08:08 Cedvega

On Sun, 2021-07-11 at 17:26 +0000, RENAUD Jean-Pierre (RDI Nancy) wrote:

What I want to do is to compute the distance to the hull of points that are outside the envelop. (these are map pixels for which I want to predict forest attributes, using a model that often have several independent variables (Xs)) As they are out of the calibration domain (out of the “convex hull” these are extrapolated according to the model that I’m using)

Brad Barber told me that I need to use 2 functions qh_findbestfacet and qh_distplane … but I don’t know how to use it in your geometry package.

It would be even more perfect if I could get not only the distance, but also the values of the intersecting points (on the Hull). I feel that not only the distance could be linked to the degree of bias, but also the direction of the vector.

@cbbarber 's message was:

Objet : Re: Help about finding distance of al point to the Hull

Many thanks for your clearly written question. What you are looking for should be available through the R package, but I am unfamiliar with their interface to Qhull. So you may need to reach to the R user groups.

You need two results for a given errant point --

  1. The facet that is closest to the point (in Qhull, qh_findbestfacet)
  2. The distance from the point to the facet (in Qhull, qh_distplane)

The distance is the distance from the point to the hyperplane for the closest facet. It is positive if the point is outside of the convex hull, as it should be in your case.

davidcsterratt avatar Aug 29 '21 13:08 davidcsterratt

The geometry package can't do what you ask at the moment, but I am now pretty confident that I can adapt the inhulln() function to do this fairly quickly.

davidcsterratt avatar Aug 29 '21 13:08 davidcsterratt

Hi David,

Something to watch for with qh_findbestfacet. It performs a directed search if the point is outside of the initial facet. Otherwise, it performs an exhaustive search.

Jean-Pierre's points are known to be outside of the convex hull, so it may be worthwhile to find a facet below the point before calling qh_findbestfacet.

I'd like to find a generic solution. CGAL has methods to reduce the set of candidate facets.

                    --Brad

At 09:01 AM 08/29/2021, David C Sterratt wrote:

The geometry package can't do what you ask at the moment, but I am now pretty confident that I can adapt the inhulln() function to do this fairly quickly.

­ You are receiving this because you were mentioned. Reply to this email directly, https://github.com/davidcsterratt/geometry/issues/59#issuecomment-907787806view it on GitHub, or https://github.com/notifications/unsubscribe-auth/ADDZ7R3JXZOZEBSLFHGXQADT7IVTTANCNFSM5CYQMXLQunsubscribe. Triage notifications on the go with GitHub Mobile for https://apps.apple.com/app/apple-store/id1477376905?ct=notification-email&mt=8&pt=524675iOS or https://play.google.com/store/apps/details?id=com.github.android&referrer=utm_campaign%3Dnotification-email%26utm_medium%3Demail%26utm_source%3DgithubAndroid.

cbbarber avatar Aug 29 '21 17:08 cbbarber

Hmm... I have been using

    facet = qh_findbestfacet(qh, point, !qh_ALL, &bestdist, &isoutside);

and then taking bestdist as the distance - this seems to work in a simple example, when the perpendicular to the plane of the facet passes through the facet itself. facet = qh_findbestfacet(qh, point, !qh_ALL, &bestdist, &isoutside);

e.g. The distance of the point (-0.5, -0.5) is found as -0.5 and (1, 1) is found as 1.412 in this simple test case: Screenshot from 2021-08-29 19-04-56 However, the point (10, -1) is found to be 6.364 from the hull; I think is this distance from the plane of the diagonal facet. I think it should be a distance of 9 from the hull, since the closest point is (1, -1). I was trying to use qh_nearvertex() to find the closest vertex, but it reports the distance as 10, and it looks like the best vertex it finds is at (0, -1) - i.e. not one of the vertices of the hull.

davidcsterratt avatar Aug 29 '21 18:08 davidcsterratt

Hi David,

For the facet closest to a point, you should call qh_findbestfacet with 'qh_ALL' instead of '!qh_ALL'. That explains your undesired result for 'the point (10, -1)'.

Extract from the documentation of qh_findbestfacet in poly2_r.c

if !bestoutside (!qh_ALL) returns first facet below point if point is inside, returns nearest, !upperdelaunay facet returns: distance to facet isoutside set if outside of facet note: If inside, qh_findbestfacet performs an exhaustive search

If you know a facet below the point (e.g., 'facetA'), then you do not need the exhaustive search branch of qh_findbestfacet. The replacement code is the first two lines of qh_findbestfacet

bestfacet= qh_findbest(qh, point, facetA, qh_ALL, !qh_ISnewfacets, qh_NOupper, bestdist, isoutside, &totpart); if (*bestdist < -qh->DISTround) { // report an error, facetA was not below the point.

                    --Brad

cbbarber avatar Aug 30 '21 22:08 cbbarber

Thanks @cbbarber. qh_ALL in qh_findbestfacet seems to work and I've pushed my work so far containing this code.

However, I've run into a conceptual problem. When a point is above a facet (i.e. the perpendicular to the hyperplane of the facet passing through the point also passes through the facet, then the distance is reported correctly. But when a point is not above a facet, e.g. the point (-2, -2) in the above example, the closest point of the hull is the nearest vertex. But how do we test if a point is above a facet in order to use the distance to the nearest vertex?

davidcsterratt avatar Sep 04 '21 20:09 davidcsterratt

The crucial changes are in this file: https://github.com/davidcsterratt/geometry/commit/5454978f943a7a33feb7e6f91b6b9b8773c6aa11#diff-18d95f23a29dae23e9245ff171c2905184d938e4f109e63072a960f8b10801e0

davidcsterratt avatar Sep 04 '21 20:09 davidcsterratt

There's a slight misunderstanding. For Qhull 'above' is relative to an interior point of the convex hull. Hyperplanes are oriented so that a point is above a facet if its signed distance to the facet's hyperplane is positive.

More properly, a point is clearly above a facet if the signed distance is greater than roundoff from the facet's outerplane. The outerplane is a positive offset from the facet's hyperplane. A point is clearly outside the convex hull if it is clearly above a facets and all of its neighbors. Nearly coplanar neighbors must be likewise tested. The definition of nearly coplanar neighbor uses empirically derived constants. After multiple merges produce wide facets, this definition may be incorrect, especially near a vertex.

A point is coplanar with a facet if it is neither clearly above nor clearly below the facet. A point is clearly below a facet, if the signed distance is less the facet's innerplane (including roundoff). The innerplane is a negative offset from the facet's hyperplane. The innerplane clearly below all of a facet's vertices.

Given the closest facet for a point, the closest vertex for the point is one of the facet's vertices or one of the vertices of a neighbor that is coplanar with point. Neighbors of a coplanar neighbor likewise need to be tested.

                                    --Brad

At 04:34 PM 09/04/2021, you wrote:

qh_ALL in qh_findbestfacet seems to work and I've pushed my work so far containing this code.

However, I've run into a conceptual problem. When a point is above a facet (i.e. the perpendicular to the hyperplane of the facet passing through the point also passes through the facet, then the distance is reported correctly. But when a point is not above a facet, e.g. the point (-2, -2) in the above example, the closest point of the hull is the nearest vertex. But how do we test if a point is above a facet in order to use the distance to the nearest vertex?

­ You are receiving this because you were mentioned. Reply to this email directly, https://github.com/davidcsterratt/geometry/issues/59#issuecomment-913036996view it on GitHub, or https://github.com/notifications/unsubscribe-auth/ADDZ7RY3DM33SKURP7S7LRDUAJ7E7ANCNFSM5CYQMXLQunsubscribe. Triage notifications on the go with GitHub Mobile for https://apps.apple.com/app/apple-store/id1477376905?ct=notification-email&mt=8&pt=524675iOS or https://play.google.com/store/apps/details?id=com.github.android&referrer=utm_campaign%3Dnotification-email%26utm_medium%3Demail%26utm_source%3DgithubAndroid.

cbbarber avatar Sep 05 '21 16:09 cbbarber