4
$\begingroup$

Suppose I have a convex hull as shown in orange. I have defined an axis for the hull by drawing an infinite line through the base, which I define, and it's centroid. I would like to know how I can find the cross sectional area of the convex hull shape at its widest point that is orthogonal to my line.

For example, in the following figures, largest cross sectional area of the orange convex hull is quite more than the purple convex hull, which is skinnier. How can I quantify this?

data = Import["filepath", {"Data"}];
cone = ConvexHullMesh[data, MeshCellStyle -> {{2, All} -> Opacity[0.2, Orange]}];
base = First[data];
centroid=RegionCentroid[cone];
axesline=Graphics3D[{Thick,InfiniteLine[{base,centroid}];
Show[cone,axesline]

thick cross sectional area

example of thinner cone

$\endgroup$
11
$\begingroup$

Here's my hack solution to the fact that RegionIntersection[ <3D mesh Region>, InfinitePlane[...]] does not evaluate. You need to get your cross section into the form of an InfinitePlane, then you find the lines which make up the intersection between the plane and the faces of the mesh region. Problem is now that these lines aren't in any order, so we take all the points making up the lines, and use FindShortestTour to put them in order, and finally make a Polygon out of these.

SeedRandom[42];
region = ConvexHullMesh[
  RandomReal[1, {14, 3}]
  ];
plane = InfinitePlane[{{1., 0.5, 1.}, {0, 0.25, 0}, {2, 0.5, 1}}];

RegionPlot3D[{region}, Axes -> True, 
  AxesLabel -> {"x", "y", "z"}]~Show~Graphics3D[{Red, plane}]

Mathematica graphics

Here we get the polygon,

crossSection = 
  RegionIntersection[plane, #] & /@ MeshPrimitives[region, 2] // 
       DeleteCases[_EmptyRegion] // 
       ReplaceAll[Line :> Sequence] // 
       Flatten[#, 1] & // 
       (#[[Last@FindShortestTour[#]]] &) // 
       Polygon;
Graphics3D[crossSection]

Mathematica graphics

Verify that this worked,

Show[
 RegionPlot3D[region, PlotStyle -> Opacity[0.5]],
 Graphics3D[{Red, crossSection}]
 ]

Mathematica graphics

Then just get the area,

Area@crossSection
(* 0.353513 *)

To get more to the exact question in the OP, if you wanted to get the cross section area when you have two points, one being the base, the other being the region centroid, you'd make up the infinite plane thusly

base = (*the origin, why not? *) {0,0,0};
point = RegionCentroid[region];
line = InfiniteLine[{base,point}];
plane = Module[{vec1,vec2,vec3,a,b,c},
    vec1=Normalize[base-point];
    vec2={a,b,c}/.First@FindInstance[{a,b,c}.vec1==0.&&Norm[{a,b,c}]==1,{a,b,c}];
    vec3=Cross[vec1,vec2];
    InfinitePlane[point,{vec2,vec3}]
];
crossSection = 
  RegionIntersection[plane, #] & /@ MeshPrimitives[region, 2] // 
       DeleteCases[_EmptyRegion] // 
       ReplaceAll[Line :> Sequence] // 
       Flatten[#, 1] & // 
       (#[[Last@FindShortestTour[#]]] &) // 
       Polygon;

Show[
    RegionPlot3D[region, PlotStyle -> Opacity[0.5]],
    Graphics3D[{
        Red, EdgeForm[Directive[Thick,Black]], crossSection,
        Yellow,Sphere[point,.02],
        Black,AbsoluteThickness[2],line}
    ],Boxed->False
 ]

Mathematica graphics

Area@crossSection
(* 0.529433 *)
$\endgroup$
  • $\begingroup$ Thank you very much! Just to clarify, I will need to define the infinite plane such that it is orthogonal to my "axis line"? $\endgroup$ – Doppler Feb 8 '17 at 19:23
  • $\begingroup$ @Doppler - see the edit, if you have 2 points, one the base and the other where you want to put the plane, then you can make the InfinitePlane pretty easily $\endgroup$ – Jason B. Feb 8 '17 at 20:27
  • $\begingroup$ @Doppler Of course there's no guarantee that the maximum cross section would contain the centroid (I think), so I may be necessary to optimize the point. $\endgroup$ – Jason B. Feb 9 '17 at 1:30
  • $\begingroup$ Hi, I tried adapting your solution but have encountered an error... "Part::pkspec1: The expression {{0.447167,0.48295,0.931799},{0.362764,0.44467,0.778682},BooleanRegion[#1&&#2&,{InfinitePlane[{{1.,0.5,1.},{0,0.25,0}... cannot be used as a part specification." Strangely enough, it was working yesterday but I have no idea what happened. Would you have any advice on this error please? $\endgroup$ – Zhao Oct 2 '17 at 10:11
  • $\begingroup$ I've checked and can confirm that your method works in Mathematica 11.0, but not 11.2, or v10 $\endgroup$ – Zhao Oct 2 '17 at 11:07

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.