Area of coplanar points in 3D - Python

Question

I got a set of coplanar points in 3D of coordinates (x,y,z) as well as the equation of the plan P they belong to. I would like to know how to get the area occupied by those points on the pane using Python. I have tried Convhull in 3d, but I get an error because all my points are indeed coplanar. Do you have any idea ?

Thanks !

I'd tackle the problem by geometrically projecting all the points onto the plane to generate the list of 2D coordinates in plane-space, then use the convex hull function in 2D only and see if that helps. — user2464424, Apr 06 '20 at 09:49
The answer to [this question](https://math.stackexchange.com/questions/3207981/caculate-area-of-polygon-in-3d) is what you are after. — DrBwts, Apr 06 '20 at 13:25

score 3 · Accepted Answer · edited Jun 20 '20 at 09:12

First I tried to reproduce your problem by creating some random, coplanar points and trying to fit a convex hull:

import numpy as np
from scipy.spatial import ConvexHull

# 4 random 3D points
points = np.random.rand(4, 3)
# A set of coplanar points where their third dimension is 1.
coplanar = np.array([np.array([x[0], x[1], 1.0]) for x in points])

hull = ConvexHull(coplanar)

This indeed generates an error:

Traceback (most recent call last):
  File "main.py", line 9, in <module>
    hull = ConvexHull(coplanar)
  File "qhull.pyx", line 2428, in scipy.spatial.qhull.ConvexHull.__init__
  File "qhull.pyx", line 357, in scipy.spatial.qhull._Qhull.__init__
scipy.spatial.qhull.QhullError: QH6154 Qhull precision error: Initial simplex is flat (facet 1 is cop
lanar with the interior point)

..... OMITTING SOME LINES ....

If the input is lower dimensional:
  - use 'QJ' to joggle the input and make it full dimensional
  - use 'Qbk:0Bk:0' to delete coordinate k from the input.  You should
    pick the coordinate with the least range.  The hull will have the
    correct topology.
  - determine the flat containing the points, rotate the points
    into a coordinate plane, and delete the other coordinates.
  - add one or more points to make the input full dimensional.

As we can see, the underlying library (qhull) gives us some suggestions in case your data is lower-dimensional. As you said, you already know that your data is coplanar and that it could be projected into a plane and be represented by 2D points.

An alternative to projecting your data, as suggested in the error message, is to use the option QJ to joggle the input (option doc). If I understood it correctly, this "joggle" makes your data "non coplanar" by introducing a random noise to your data, allowing the optimization to proceed.

hull = ConvexHull(coplanar, qhull_options="QJ")
hull.area
> 0.3100618849870332

In summary:

find the angle between the plane and the unit xy plane, and rotate your data if you want a precise answer
use the option "QJ" if you want a quick solution and high precision is not a problem (I actually don't know how far away the answer can go from the truth)

Update

From Qhull FAQ:

The area is the area of the surface of the convex hull, while the volume is the total volume of the convex hull.

In 2-d, the convex hull is a polygon. Its surface is the edges of a polygon. So in 2-d, the 'area' is the length of the polygon's edges, while the 'volume' is the area of the polygon.

From these statements, I interpreted that you should use:

hull.area / 2 for 3D points (your original points with the option QJ)
hull.volume for 2D points (if you rotate your points and get rid of one dimension)

To help clarify the usage of ConvexHull, I used a simpler example

square2D = np.array([
    [0.0, 0.0],
    [2.0, 0.0],
    [2.0, 2.0],
    [0.0, 2.0]
])
hull = ConvexHull(square2D)
hull.area
> 8.0
hull.volume
> 4.0

The results are consistent with the documentation.

Now, to understand the effects of the option QJ, I just added one dimension to the previous points:

coplanar = np.array([
    [0.0, 0.0, 0.0],
    [2.0, 0.0, 0.0],
    [2.0, 2.0, 0.0],
    [0.0, 2.0, 0.0]
])
hull = ConvexHull(coplanar, qhull_options="QJ")
hull.area
> 8.000000000306578
hull.volume
> 5.719619827513867e-11

For this last example, I imagine the convex hull as an approximately "flat" 3D object (because of the introduced noise) with two big facets. This is why I thought that you could use hull.area / 2 in this case.

Thank you for your answer ! It indeed is my problem, I will check on what QJ is then ! For the prjection, would you know any python module able to make a projection ? (what I would need I I get convhull and the error right would be local coordinates (u,v) from global coordinates of my coplanar points (x,y,z)) — cazator, Apr 06 '20 at 12:04
@cazator, this is a quick thought, I may need to elaborate it better (my math is really rusty). First, you need to find the angle between the plane and the unit coordinate (you may use the tangent equation); second, you could use `scipy.spatial.transform.Rotation` to rotate the point vectors in order to have `z = 0`; finally, you can reshape your final point vectors and remove the third dimension. — boechat107, Apr 06 '20 at 12:34
be careful with projection as the projected polyon will almost certainly have a different area than the original. A better approach would be to calculate the rotation matrix that would tranform your plane to the xy plane. — DrBwts, Apr 06 '20 at 13:08
Hi thank you for your answers. I went into the math, and computed the local coordinates (u,v) of my points on the plane, before doing ConvHull on the (u,v) coordinates. I get the following points and ConvHull results : https://imgur.com/a/BhPIBVG. I checked, the coordinates seem to be OK in terms of orders of magnitude and pattern. ConvHull seems to work well for the shape also. But the area given by ConvHull(.).area is around 4000, which seems small as I expected sthg 2 orders of magnitude higher. Do you guys have any idea on why the area is so small ? Even if the algorithm seems to work ? — cazator, Apr 06 '20 at 14:30
@DrBwts is totally correct. I've used the wrong term trying to say the same thing. I'll edit the answer. — boechat107, Apr 06 '20 at 15:05
@cazator, what about the area calculated with the option `QJ`? It may give you an approximate answer or a reference value. Besides that, it might be a problem with the rotation. — boechat107, Apr 06 '20 at 15:10
I am going to try with QJ, it sounds interesting with the link you provided. I don't think it is a problem with the rotation itself. Indeed if we consider the rotated set of points as the original one, the area should still be higher than the one returned. I found there (https://stackoverflow.com/questions/52403793/perimeter-of-a-2d-convex-hull-in-python) that for the area one should actually type Convhull(...).volume and not .area when in 2D, but I don't get why. And with .volume, the area returned is of the good order of magnitude (assuming the points on the picture are in a circle ) — cazator, Apr 06 '20 at 15:36
@cazator, you are right, there is a difference between using for 2D and 3D points. I've expanded the answer to include this information. — boechat107, Apr 07 '20 at 10:51

Area of coplanar points in 3D - Python

1 Answers1

Update