Polygon envelope of a 3D array

I've written the attached function to find the polygon envelope of structured data in a 3D array. It works, but is very slow for large data sets. Does anyone know of a way to make it substantially faster?

Also, it is unfortunately limited to 3D array input (the Vert/Faces) approach enforces this limit. Does anyone see a way to generalize this to ND data?

If this were fast and generalized, I think it would make a nice addition to the File Exchange...

[Andrew Newell]

The function is not attached.

[Rob McDonald]

Sorry. It is now.

Answer 1

[Image Analyst]

Any reason why you didn't just use the built in convhulln()?

[Rob McDonald]

The shape may not be convex.

[Image Analyst]

Then, can you then give us a high level description of how you're defining "envelope"? If your envelope can "dip inside" the convex hull and attach to arbitrary vertices inside, then there are tons of vertices it could connect to so how do you decide which is the "right" one? Why aren't they all right? It seems you could have hundreds or thousands of possible envelopes.

[Rob McDonald]

The envelope is implied by the structure of the points in the array. They have been created by evaluating continuous functions of three variables over a certain domain.

The envelope approximately contains all possible solutions to the functions across the domain. The approximation improves as the number of sample points across the domain increases.

In the sample code, this is accomplished by treating neighboring points in the 3D arrays as topological cubes. Each of the faces of the cubes are extracted as polygons projections. The envelope is cumulatively constructed by taking the Boolean union of the polygons with the envelope.

In the below example (from a study I did last year), the white envelope represents the complete torque/rpm envelope required to turn a fixed pitch propeller through a wide range of throttle setting (also rpm), altitude, and airspeed. The peak represents sea-level, static, full-throttle operation. The curve back to the origin corresponds to decreasing throttle at sea-level static conditions. The jagged bottom is caused by NaN's returned by the propeller model along its stall boundary.

Similarly (though created from 2D arrays), the grey envelope represents the complete torque/rpm envelope for that propeller with the throttle adjusted to achieve equilibrium (Thrust=Drag) flight for an assumed aircraft (across all possible airspeeds and altitudes where T=D can be achieved for this aircraft/prop pair). In this case, the curved left boundary of the envelope is formed by an airspeed sweep at low altitude. As airspeed increases, the torque decreases and then increases while the speed increases monotonically.

These envelopes may be used to consider the size and performance of a motor or engine to turn the prop.

The most likely way to speed this up would be to find a way to vectorize the call to polybool to combine more than one poly at a time.

The documentation says you can pass a nan-delimited list (or cell array) of polys, but I'm not entirely sure of the semantics of (A U B) when A or B itself a list of polys.

How to group the input polys is also a question -- all at once, one cube's worth at a time, a strip through the data, a slice, etc?

I currently add all six faces of each cube -- so I end up adding every interior face twice. Slightly more clever looping would get about 2x speedup.

More exotically, perhaps a re-entrant version of polybool that keeps the envelope in memory as polygons are added. Perhaps a spatial data structure to pre-check the polys before the call to polybool.

I was tempted to approach this by focusing on the boundaries of the domain, but there are cases (like the bottom edge of the grey envelope) where the envelope is formed by the interior of the domain.

[Rob McDonald]

Here is an updated version that avoids duplicate checks on interior faces. As expected, this achieves about a 2x speedup on a large problem.

Passing a list of polys doesn't automatically self-union the list, so any approach using lists would have to make sure they overlap. A sequence of i-strips then j-strips then k-strips may work.

However, figuring out when/how to start/stop each strip if NaN is encountered in the data is not obvious.