The Perl Toolchain Summit needs more sponsors. If your company depends on Perl, please support this very important event.

NAME

Math::Geometry::Delaunay - Quality Mesh Generator and Delaunay Triangulator

VERSION

Version 0.21

SYNOPSIS



input vertices

Delaunay triangulation

Voronoi diagram
    use Math::Geometry::Delaunay qw(TRI_CCDT);

    # generate Delaunay triangulation
    # and Voronoi diagram for a point set

    my $point_set = [ [1,1], [7,1], [7,3],
                      [3,3], [3,5], [1,5] ];

    my $tri = new Math::Geometry::Delaunay();
    $tri->addPoints($point_set);
    $tri->doEdges(1);
    $tri->doVoronoi(1);
    
    # called in void context
    $tri->triangulate();
    # populates the following lists

    $tri->elements(); # triangles
    $tri->nodes();    # points
    $tri->edges();    # triangle edges
    $tri->vnodes();   # Voronoi diagram points
    $tri->vedges();   # Voronoi edges and rays


input PSLG

output mesh

something interesting
extracted from topology
    # quality mesh of a planar straight line graph
    # with cross-referenced topological output

    my $tri = new Math::Geometry::Delaunay();
    $tri->addPolygon($point_set);
    $tri->minimum_angle(23);
    $tri->doEdges(1);

    # called in scalar context
    my $mesh_topology = $tri->triangulate(TRI_CCDT);
    # returns cross-referenced topology

    # make two lists of triangles that touch boundary segments

    my @tris_with_boundary_segment;
    my @tris_with_boundary_point;

    foreach my $triangle (@{$mesh_topology->{elements}}) {
            my $nodes_on_boundary_count = ( 
                grep $_->{marker} == 1,
                @{$triangle->{nodes}} 
                );
            if ($nodes_on_boundary_count == 2) {
                push @tris_with_boundary_segment, $triangle;
                }
            elsif ($nodes_on_boundary_count == 1) {
                push @tris_with_boundary_point, $triangle;
                }
            }
            

DESCRIPTION

This is a Perl interface to the Jonathan Shewchuk's Triangle library.

"Triangle generates exact Delaunay triangulations, constrained Delaunay triangulations, conforming Delaunay triangulations, Voronoi diagrams, and high-quality triangular meshes. The latter can be generated with no small or large angles, and are thus suitable for finite element analysis." -- from http://www.cs.cmu.edu/~quake/triangle.html

EXPORTS

Triangle has several option switches that can be used in different combinations to choose a class of triangulation and then configure options within that class. To clarify the composition of option strings, or just to give you a head start, a few constants are supplied to configure different classes of mesh output.

    TRI_CONSTRAINED  = 'Y'    for "Constrained Delaunay"
    TRI_CONFORMING   = 'Dq0'  for "Conforming Delaunay"
    TRI_CCDT         = 'q'    for "Constrained Conforming Delaunay"
    TRI_VORONOI      = 'v'    to generate the Voronoi diagram

For an illustration of these terms, see: http://www.cs.cmu.edu/~quake/triangle.defs.html

CONSTRUCTOR

new

The constructor returns a Math::Geometry::Delaunay object.

    my $tri = Math::Geometry::Delaunay->new();

MESH GENERATION

triangulate

Run the triangulation with specified options, and either populate the object's output lists, or return a hash reference giving access to a cross-referenced representation of the mesh topology.

Common options can be set prior to calling triangulate. The full range of Triangle's options can also be passed to triangulate as a string, or list of strings. For example:

    my $tri = Math::Geometry::Delaunay->new('pzq0eQ');

    my $tri = Math::Geometry::Delaunay->new(TRI_CCDT, 'q15', 'a3.5');

Triangle's command line switches are documented here: http://www.cs.cmu.edu/~quake/triangle.switch.html

list output

After triangulate is invoked in void context, the output mesh data can be retrieved from the following methods, all of which return a reference to an array.

    $tri->triangulate(); # void context - no return value requested
    # output lists now available
    $points  = $tri->nodes();    # array of vertices
    $tris    = $tri->elements(); # array of triangles
    $edges   = $tri->edges();    # all the triangle edges
    $segs    = $tri->segments(); # the PSLG segments
    $vpoints = $tri->vnodes();   # points in the voronoi diagram
    $vedges  = $tri->vedges();   # edges in the voronoi diagram

Data may not be available for all lists, depending on which option switches were used. By default, nodes and elements are generated, while edges are not.

The members of the lists returned have these formats:

    nodes:    [x, y, < zero or more attributes >, < boundary marker >]

    elements: [[x0, y0], [x1, y1], [x2, y2],
                < another three vertices, if "o2" switch used >,
                < zero or more attributes >
                ]
    edges:    [[x0, y0], [x1, y1], < boundary marker >]

    segments: [[x0, y0], [x1, y1], < boundary marker >]

    vnodes:   [x, y, < zero or more attributes >]

    vedges:   [< vertex or vector >, < vertex or vector >, < ray flag >]

Boundary markers are 1 or 0. An edge or segment with only one end on a boundary has boundary marker 0.

The ray flag is 0 if the edge is not a ray, or 1 or 2, to indicate which vertex is actually a unit vector indicating the direction of the ray.

Import of the mesh data from the C data structures will be deferred until actually requested from the list fetching methods above. For speed and lower memory footprint, access only what you need, and consider suppressing output you don't need with option switches.

topological output

When triangulate is invoked in scalar or array context, it returns a hash ref containing the cross-referenced nodes, elements, edges, and PSLG segments of the triangulation. In array context, with the "v" switch enabled, the Voronoi topology is the second item returned.

    my $topology = $tri->triangulate();

    $topology now looks like this:
    
    {
    nodes    => [
                  { # a node
                  point      => [x0, x1],
                  edges      => [edgeref, ...],
                  segments   => [edgeref, ...], # a subset of edges
                  elements   => [elementref, ...],
                  marker     => 1 or 0 or undefined, # boundary marker
                  attributes => [attr0, ...]
                  },
                  ... more nodes like that
                    
                ],
    elements => [
                  { # a triangle
                  nodes      => [noderef0, noderef1, noderef2],
                  edges      => [edgeref0, edgeref1, edgeref2],
                  neighbors  => [neighref0, neighref1, neighref2],
                  attributes => [attrib0, ...]
                  },
                  ... more triangles like that
                ],
    edges    => [
                  {
                  nodes    => [noderef0, noderef1], # only one for a ray
                  elements => [elemref0, elemref1], # one if on boundary
                  vector   => undefined or [x, y],  # ray direction
                  marker   => 1 or 0 or undefined,  # boundary marker
                  index    => <integer> # edge's index in edge list
                  },
                  ... more edges like that
                    
                ],
    segments => [
                  {
                  nodes    => [noderef0, noderef1],
                  elements => [elemref0, elemref1], # one if on boundary
                  marker   => 1 or 0 or undefined   # boundary marker
                  },
                  ... more segments
                ]
    }

cross-referencing Delaunay and Voronoi

Corresponding Delaunay triangles and Voronoi nodes have the same index number in their respective lists.

In the topological output, any element in a triangulation has a record of its own index number that can by used to look up the corresponding node in the Voronoi diagram topology, or vice versa, like so:

    ($topo, $voronoi_topo) = $tri->triangulate('v');

    # get a triangle reference where the index is not obvious
    
    $element = $topo->{nodes}->[-1]->{elements}->[-1];
    
    # this gets a reference to the corresponding node in the Voronoi diagram
    
    $voronoi_node = $voronoi_topo->{nodes}->[$element->{index}];

Corresponding edges in the Delaunay and Voronoi outputs have the same index number in their respective edge lists.

In the topological output, any edge in a triangulation has a record of its own index number that can by used to look up the corresponding edge in the Voronoi diagram topology, or vice versa, like so:

    ($topo, $voronoi_topo) = $tri->triangulate('ev');
    
    # get an edge reference where it's not obvious what the edge's index is
    
    $delaunay_edge = $topo->{nodes}->[-1]->{edges}->[-1];
    
    # this gets a reference to the corresponding edge in the Voronoi diagram
    
    $voronoi_edge = $voronoi_topo->{edges}->[$delaunay_edge->{index}];

METHODS TO SET SOME Triangle OPTIONS

area_constraint

Corresponds to the "a" switch.

With one argument, sets the maximum triangle area constraint for the triangulation. Returns the value supplied. With no argument, returns the current area constraint.

Passing -1 to area_constraint() will disable the global area constraint.

minimum_angle

Corresponds to the "q" switch.

With one argument, sets the minimum angle allowed for triangles added in the triangulation. Returns the value supplied. With no argument, returns the current minimum angle constraint.

Passing -1 to minimum_angle() will cause the "q" switch to be omitted from the option string.

doEdges, doVoronoi, doNeighbors

These methods simply add or remove the corresponding letters from the option string. Pass in a true or false value to enable or disable. Invoke with no argument to read the current state.

quiet, verbose

Triangle prints a basic summary of the meshing operation to STDOUT unless the "Q" switch is present. This module includes the "Q" switch by default, but you can override this by passing a false value to quiet().

If you would like to see even more output regarding the triangulation process, there are are three levels of verbosity configurable with repeated "V" switches. Passing a number from 1 to 3 to the verbose() method will enable the corresponding level of verbosity.

METHODS TO ADD VERTICES AND SEGMENTS

addVertices, addPoints

Takes a reference to an array of vertices, each vertex itself an reference to an array containing two coordinates and zero or more attributes. Attributes are floating point numbers.

    # vertex format
    # [x, y, < zero or more attributes as floating point numbers >]

    $tri->addPoints([[$x0, $y0], [$x1, $y1], ... ]);

Use addVertices to add vertices that are not part of a PSLG. Use addPoints to add points that are not part of a polygon or polyline. In other words, they do the same thing.

addSegments

Takes a reference to an array of segments.

    # segment format
    # [[$x0, $y0], [$x1, $y1]]

    $tri->addSegments([ $segment0, $segment1, ... ]);

If your segments are contiguous, it's better to use addPolyline, or addPolygon.

This method is provided because some point and polygon processing algorithms result in segments that represent polygons, but list the segments in a non-contiguous order, and have shared vertices repeated in each segment's record.

The segments added with this method will be checked for duplicate vertices, and references to these will be merged.

Triangle can handle duplicate vertices, but we would rather not feed them in on purpose.

addPolyline

Takes a reference to an array of vertices describing a curve. Creates PSLG segments for each pair of adjacent vertices. Adds the new segments and vertices to the triangulation input.

    $tri->addPolyline([$vertex0, $vertex1, $vertex2, ...]);

addPolygon

Takes a reference to an array of vertices describing a polygon. Creates PSLG segments for each pair of adjacent vertices and creates and additional segment linking the last vertex to the first,to close the polygon. Adds the new segments and vertices to the triangulation input.

    $tri->addPolygon([$vertex0, $vertex1, $vertex2, ...]);

addHole

Like addPolygon, but describing a hole or concavity - an area of the output mesh that should not be triangulated.

There are two ways to specify a hole. Either provide a list of vertices, like for addPolygon, or provide a single vertex that lies inside of a polygon, to identify that polygon as a hole.

    # first way
    $tri->addHole([$vertex0, $vertex1, $vertex2, ...]);

    # second way
    $tri->addPolygon( [ [0,0], [1,0], [1,1], [0,1] ] );
    $tri->addHole( [0.5,0.5] );

Hole marker points can also be used, in combination with the "c" option, to cause or preserve concavities in a boundary when Triangle would otherwise enclose a PSLG in a convex hull.

addRegion

Takes a polygon describing a region, and an attribute or area constraint. With both the "A" and "a" switches in effect, three arguments allow you to specify both an attribute and an optional area constraint.

The first argument may alternately be a single vertex that lies inside of another polygon, to identify that polygon as a region.

To be used in conjunction with the "A" and "a" switches.

    # with the "A" switch
    $tri->addRegion(\@polygon, < attribute > );
    
    # with the "a" switch
    $tri->addRegion(\@polygon, < area constraint > );

    # with both "Aa"
    $tri->addRegion(\@polygon, < attribute >, < area constraint > );

If the "A" switch is used, each triangle generated within the bounds of a region will have that region's attribute added to the end of the triangle's attributes list, while each triangle not within a region will have a "0" added to the end of its attribute list.

If the "a" switch is used without a number following, each triangle generated within the bounds of a region will be subject to that region's area constraint.

If the "A" or "a" switches are not in effect, addRegion has the same effect as addPolygon.

METHODS TO ACCESS OUTPUT LISTS

The following methods retrieve the output lists after the triangulate method has been invoked in void context.

Triangle's output data is not imported from C to Perl until one of these methods is invoked, and then only what's needed to construct the list requested. So there may be a speed or memory advantage to accessing the output in this way - only what you need, when you need it.

The methods prefixed with "v" access the Voronoi diagram nodes and edges, if one was generated.

nodes

Returns a reference to a list of nodes (vertices or points).

    my $pointlist = $tri->nodes();    # retrieve nodes/vertices/points
    

The nodes in the list have this structure:

    [x, y, < zero or more attributes >, < boundary marker >]

elements

Returns a reference to a list of elements.

    $triangles  = $tri->elements(); # retrieve triangle list

The elements in the list have this structure:

    [[x0, y0], [x1, y1], [x2, y2],
     < another three vertices, if "o2" switch used >
     < zero or more attributes >
    ]

segments

Returns a reference to a list of segments.

    $segs  = $tri->segments(); # retrieve the PSLG segments

The segments in the list have this structure:

    [[x0, y0], [x1, y1], < boundary marker >]

edges

Returns a reference to a list of edges.

    $edges  = $tri->edges();    # retrieve all the triangle edges

The edges in the list have this structure:

    [[x0, y0], [x1, y1], < boundary marker >]

Note that the edge list is not produced by default. Request that it be generated by invoking doEdges(1), or passing the 'e' switch to triangulate().

vnodes

Returns a reference to a list of nodes in the Voronoi diagram.

    $vpointlist = $tri->vnodes();   # retrieve Voronoi vertices

The Voronoi diagram nodes in the list have this structure:

    [x, y, < zero or more attributes >]

vedges

Returns a reference to a list of edges in the Voronoi diagram. Some of these edges are actually rays.

    $vedges = $tri->vedges();   # retrieve Voronoi diagram edges and rays 

The Voronoi diagram edges in the list have this structure:

    [< vertex or vector >, < vertex or vector >, < ray flag >]

If the edge is a true edge, the ray flag will be 0. If the edge is actually a ray, the ray flag will either be 1 or 2, to indicate whether the the first, or second vertex should be interpreted as a direction vector for the ray.

UTILITY FUNCTIONS

to_svg

This function is meant as a development and debugging aid, to "dump" the geometric data structures specific to this package to a graphical representation. Takes key-value pairs to specify topology hashes, output file, image dimensions, and styles for the elements in the various output lists.

The topology hash input for the topo or vtopo keys is just the hash returned by triangulate. The value for the file key is a file name string. Omit file to print to STDOUT. For size, provide and array ref with width and height, in pixels. For output list styles, keys correspond to the output list names, and values consist of references to arrays containing style configurations, as demonstrated below.

Only geometry that has a style configuration will be displayed. The following example includes everything. To display a subset, just omit any of the style configuration key-value pairs.

    ($topo, $vtopo) = $tri->triangulate('ve');

    to_svg( topo  => $topo,
            vtopo => $vtopo,
            
            file => "enchilada.svg",    # omit for STDOUT
            size => [800, 600],         # width, height in pixels
            
            #                     line width or   optional
            #         svg color   point radius    extra CSS
            
            nodes    => ['black'  ,   0.3],
            edges    => ['#CCCCCC',   0.7],
            segments => ['blue'   ,   0.9,     'stroke-dasharray:1 1;'],
            elements => ['pink']  , # string or callback; see below

            # these require Voronoi input (vtopo)

            vnodes   => ['purple' ,   0.3],
            vedges   => ['#FF0000',   0.7],
            vrays    => ['purple' ,   0.6],
            circles  => ['orange' ,   0.6],
            
          );

Note that for display purposes vedges does not include the infinite rays in the Voronoi diagram. To see the complete Voronoi diagram, including segments representing the infinite rays, you should include style configuration for the vrays key, as in the example above.

Elements (triangles) only need one style config entry, for color. (An optional second entry would be a string for additional CSS.) In this case, the first entry can also be a reference to a callback function. A reference to the triangle being processed for display will be passed to the callback function. Therefore the callback function can determine a color based on any features or relationships of that triangle.

Typically you might color each triangle according to the region it's in, by using Triangle's 'A' switch, and then reading the region attribute from the last item in the triangle's attribute list.

    my $region_colors_callback = sub {
        my $tri_ref = shift;
        return ('gray','blue','green')[$tri_ref->{attributes}->[-1]];
        };

But any other data accessible through the triangle reference can be used to calculate a color. For instance, the triangle's three nodes can carry any number of attributes, which are interpolated during mesh generation. You might shade each triangle according to the average of a node attribute.

    my $tri_nodes_average_callback = sub {
        my $tri_ref = shift;
        my $sum = 0;
        # calculate average of the eighth attribute in all nodes
        foreach my $node (@{$tri_ref->{nodes}}) {
            $sum += $node->{attributes}->[7];
            }
        return &attrib_val_to_grayscale_hexcode( $sum / 3 );
        };

mic_adjust


Voronoi edges (blue)
as a poor medial
axis approximation

improved approximation
after calling mic_adust()

Warning: not yet thoroughly tested; may move elsewhere

One use of the Voronoi diagram of a tessellated polygon is to derive an approximation of the polygon's medial axis by pruning infinite rays and perhaps trimming or refining remaining branches. The approximation improves as intervals between sample points on the polygon become shorter. But it's not always desirable to multiply the number of polygon points to achieve short intervals.

At any point on the true medial axis, there is a maximum inscribed circle, with it's center on the medial axis, and tangent to the polygon in at least two places.

The mic_adjust() function moves each Voronoi node so that it becomes the center of a circle that is tangent to the polygon at two points. In simple cases this is a maximum inscribed circle, and the point is on the medial axis. And when it's not, it still should be a much better approximation than the original point location. The radius to the tangent on the polygon is stored with the updated Voronoi node.

After calling mic_adjust(), the modified Voronoi topology can be used as a list of maximum inscribed circles, from which can be derive a straighter, better medial axis approximation, without having to increase the number of sample points on the polygon.

    ($topo, $voronoi_topo) = $tri->triangulate('e');

    mic_adjust($topo, $voronoi_topo); # modifies $voronoi_topo in place
    
    foreach my $node (@{$voronoi_topo->{nodes}}) {
        $mic_center = $node->{point};
        $mic_radius = $node->{radius};
        ...
        }

Constructing a true medial axis is much more involved - a subject for a different module. Until that module appears, running topology through mic_adjust() and then walking and pruning the Voronoi topology might help fill the gap.

API STATUS

Currently Triangle's option strings are exposed to give more complete access to its features. More of these options, and perhaps certain common combinations of them, will likely be wrapped in method-call getter-setters. I would prefer to preserve the ability to use the option strings directly, but it may be better at some point to hide them completely behind a more descriptive interface.

AUTHOR

Michael E. Sheldrake, <sheldrake at cpan.org>

Triangle's author is Jonathan Richard Shewchuk

BUGS

Please report any bugs or feature requests to

bug-math-geometry-delaunay at rt.cpan.org

or through the web interface at

http://rt.cpan.org/NoAuth/ReportBug.html?Queue=Math-Geometry-Delaunay

I will be notified, and then you'll automatically be notified of progress on your bug as I make changes.

SUPPORT

You can find documentation for this module with the perldoc command.

    perldoc Math::Geometry::Delaunay

You can also look for information at:

ACKNOWLEDGEMENTS

Thanks go to Far Leaves Tea in Berkeley for providing oolongs and refuge, and a place for paths to intersect.

LICENSE AND COPYRIGHT

Copyright 2013 Micheal E. Sheldrake.

This Perl binding to Triangle is free software; you can redistribute it and/or modify it under the terms of either: the GNU General Public License as published by the Free Software Foundation; or the Artistic License.

See http://dev.perl.org/licenses/ for more information.

Triangle license

Triangle by Jonathan Richard Shewchuk, copyright 2005, includes the following notice in the C source code. Please refer to the C source, included in with this Perl module distribution, for the full notice.

    This program may be freely redistributed under the condition that the
    copyright notices (including this entire header and the copyright
    notice printed when the `-h' switch is selected) are not removed, and
    no compensation is received.  Private, research, and institutional
    use is free.  You may distribute modified versions of this code UNDER
    THE CONDITION THAT THIS CODE AND ANY MODIFICATIONS MADE TO IT IN THE
    SAME FILE REMAIN UNDER COPYRIGHT OF THE ORIGINAL AUTHOR, BOTH SOURCE
    AND OBJECT CODE ARE MADE FREELY AVAILABLE WITHOUT CHARGE, AND CLEAR
    NOTICE IS GIVEN OF THE MODIFICATIONS.  Distribution of this code as
    part of a commercial system is permissible ONLY BY DIRECT ARRANGEMENT
    WITH THE AUTHOR.  (If you are not directly supplying this code to a
    customer, and you are instead telling them how they can obtain it for
    free, then you are not required to make any arrangement with me.)