# How to calculate volume of convex hull and volume of a 3D object

I have a set of random 3D data points.

How can I calculate the volume of the convex hull?

-
What is your question? –  R Hall Jun 16 '12 at 21:53
The "area" of fig01 is not well defined. The usual definiton of "Area" applies to 2D surfaces, not to point sets –  belisarius Jun 16 '12 at 21:58
For your new question, see this. –  Ｊ. Ｍ. Jun 17 '12 at 4:29

This following code that uses TetGen we will compute the volume of the convex hull.

Needs["TetGenLink"];
TetraMaker[pts_, surface_, TetGenString_?StringQ] :=
Module[{inInst, outInst, coords, surface1, meshElements, facets},
inInst = TetGenCreate[];
TetGenSetPoints[inInst, pts];
facets = Partition[surface, 1];
TetGenSetFacets[inInst, facets];
outInst = TetGenTetrahedralize[inInst, TetGenString];
coords = TetGenGetPoints[outInst];
surface1 = TetGenGetFaces[outInst];
meshElements = TetGenGetElements[outInst];
{coords, surface1, meshElements}
];
TetrahedraVolume = Compile[{{coords, _Real, 2}, {elements, _Integer, 1}},
Block[{p},
p = coords[[elements]];
1/6*Abs[Det[p[[ {1, 2, 3}]] - p[[{2, 3, 4}]]]]
], RuntimeAttributes -> {Listable}, RuntimeOptions -> "Speed"
];


Some 3D data

data3D = RandomReal[{0, 10}, {65, 3}];


Compute the convex hull and then call the above function to form the tetrahedralization. Then call TetrahedraVolume to compute the volume.

{pts, surface} = TetGenConvexHull[data3D];
{coords, surface1, meshElements} = TetraMaker[pts, surface, "pqa2.8"];
Total[ TetrahedraVolume[coords, meshElements]]


505.135

You can use this to compute surface area of a triangulated 3D geometry

TriangleArea[pts_List?(Length[#] == 3 &)] :=
Norm[Cross[pts[[2]] - pts[[1]], pts[[3]] - pts[[1]]]]/2
TriangleArea[{pts[[#[[1]]]], pts[[#[[2]]]], pts[[#[[3]]]]}] & /@ surface // Total


337.121

I compute the area of the triangles separately and adding them gives me the area of the surface that defines the convex hull.

In the above picture first you see the convex hull in black lines. The middle one shows the blue surface mesh created by TetGen during tetrahedralization. In the last one you can see the cell volumes of the tetrahedrons that discretize the volume of the convex hull in different random colors. We get the total volume by adding the volumes of these tetrahedrons.

## Volume for Convex Hull of the point cloud in your data

data3D = Import["http://dl.dropbox.com/u/68983831/object.vtk", "VertexData"];
{pts, surface} = TetGenConvexHull[data3D];
{coords, surface1, meshElements} = TetraMaker[pts, surface, "pqa2.8"];
Total[ TetrahedraVolume[coords, meshElements]]


3120.05

## Volume for 3D Geometry in your data

surface = Import["http://dl.dropbox.com/u/68983831/object.vtk", "PolygonData"];
{coords, surface1, meshElements} = TetraMaker[pts, surface, "pqa.8"];
Total[TetrahedraVolume[coords, meshElements]]


1622.23

Code for Graphics

p1 = Graphics3D[{{Red, PointSize[0.03], Point[data3D]}, {Yellow,
Opacity[.8], EdgeForm[{Thick, Black}],
GraphicsComplex[pts, Polygon[surface]]}}, Boxed -> False, Axes -> True];
p2 = Graphics3D[{{Red, PointSize[0.03],
Point[data3D]}, {EdgeForm[{Thick, Black}], FaceForm[None],
GraphicsComplex[pts, Polygon[surface]]}, {Yellow, Opacity[.8],
EdgeForm[Blue], GraphicsComplex[coords, Polygon[surface1]]}},
Boxed -> False, Axes -> True];
p3 = Graphics3D[{{Red, PointSize[0.03], Point[data3D]},
Table[With[{p =
RGBColor[RandomReal[{0, 1}, 3]]}, {Blend[{Lighter[Yellow, .8],
Red, Yellow}, i/Length@meshElements], Glow[Darker[p, .5]],
Specularity[White, 20], Opacity[.2], EdgeForm[p],
GraphicsComplex[coords,
Polygon@Partition[meshElements[[i]], 3, 1, 1]]}], {i, 1,
Length@meshElements}]}, Boxed -> False, Axes -> True];
GraphicsGrid[Transpose@{{p1, p2, p3}}, Spacings -> {0, 0},ImageSize -> 600]

-
There's a less messy way to compute surface area: Total[Flatten[Normal[GraphicsComplex[pts, Polygon[surface]]]] /. Polygon -> TriangleArea] –  Ｊ. Ｍ. Jun 17 '12 at 3:54

Taken together, PlatoManiac's answer and J.M.'s comment show how to do this. Here's code:

Needs["TetGenLink"];
data3D = Import["http://dl.dropbox.com/u/68983831/object.vtk", "VertexData"];
{pts, surface} = TetGenConvexHull[data3D];
volumeContribution[{v1_, v2_, v3_}] := With[
{n = Cross[v2 - v1, v3 - v1]}, Abs[v1.n]/2];
triangles = pts[[#]] & /@ surface;
Total[volumeContribution /@ triangles]/3


20491.

-
I think that should be v1.n on the fifth line, not Abs[v1.n]. –  Heike Jun 17 '12 at 19:10
@ MarkMcClure: Thanks +1 –  mathew Jun 18 '12 at 3:32

Use RegionMeasure and ConvexHullMesh. Example:

points = RandomReal[{-1, 1}, {10, 3}]

RegionMeasure@ConvexHullMesh[points]

(* ==> 0.948199 *)

-
This is clearly a forgotten gem +1. –  bobthechemist Sep 23 at 18:40