Triangulation Statistics
DelaunayTriangulation.IndividualTriangleStatistics
— TypeIndividualTriangleStatistics{T}
Struct storing statistics of a single triangle.
Fields
area::T
: The area of the triangle.lengths::NTuple{3,T}
: The lengths of the edges of the triangle, given in sorted order.circumcenter::NTuple{2,T}
: The circumcenter of the triangle.circumradius::T
: The circumradius of the triangle.angles::NTuple{3, T}
: The angles of the triangle, given in sorted order.radius_edge_ratio::T
: The ratio of the circumradius to the shortest edge length.edge_midpoints::NTuple{3,NTuple{2,T}}
: The midpoints of the edges of the triangle.aspect_ratio::T
: The ratio of the inradius to the circumradius.inradius::T
: The inradius of the triangle.perimeter::T
: The perimeter of the triangle.centroid::NTuple{2,T}
: The centroid of the triangle.offcenter::NTuple{2,T}
: The offcenter of the triangle with radius-edge ratio cutoffβ=1
. See this paper.sink::NTuple{2,T}
: The sink of the triangle relative to the parent triangulation. See this paper.
Constructors
The constructor is
IndividualTriangleStatistics(p, q, r, sink = (NaN, NaN))
where p
, q
, and r
are the coordinates of the triangle given in counter-clockwise order. sink
is the triangle's sink. This must be provided separately since it is only computed relative to a triangulation, and so requires vertices rather than coordinates; see triangle_sink
.
Extended help
The relevant functions used for computing these statistics are
DelaunayTriangulation.triangle_area
— Functiontriangle_area(p, q, r) -> Number
Returns the signed area of a triangle (p, q, r)
. The area is positive if (p, q, r)
is positively oriented.
triangle_area(ℓ₁²::Number, ℓ₂²::Number, ℓ₃²::Number) -> Number
Compute the area of a triangle given the squares of its edge lengths. The edges should be sorted so that ℓ₁² ≤ ℓ₂² ≤ ℓ₃²
. If there are precision issues that cause the area to be negative, then the area is set to zero.
DelaunayTriangulation.triangle_circumradius
— Functiontriangle_circumradius(A, ℓmin², ℓmed², ℓmax²) -> Number
Computes the circumradius of a triangle with area A
and squared edge lengths ℓmin² ≤ ℓmed² ≤ ℓmax²
. The circumradius is given by
\[r = \dfrac{\ell_{\min}\ell_{\text{med}}\ell_{\max}}{4A}.\]
triangle_circumradius(p, q, r) -> Number
Computes the circumradius of the triangle with coordinates (p, q, r)
.
DelaunayTriangulation.triangle_perimeter
— Functiontriangle_perimeter(ℓmin::Number, ℓmed::Number, ℓmax::Number) -> Number
Computes the perimeter of a triangle with edge lengths ℓmin ≤ ℓmed ≤ ℓmax
. The perimeter is given by
\[P = \ell_{\min} + \ell_{\text{med}} + \ell_{\max}.\]
triangle_perimeter(p, q, r) -> Number
Computes the perimeter of the triangle with coordinates (p, q, r)
.
DelaunayTriangulation.triangle_inradius
— Functiontriangle_inradius(A, perimeter) -> Number
Computes the inradius of a triangle with area A
and perimeter perimeter
. The inradius is given by
\[r_i = \dfrac{2A}{P},\]
where $P$ is the perimeter
.
triangle_inradius(p, q, r) -> Number
Computes the inradius of the triangle with coordinates (p, q, r)
.
DelaunayTriangulation.triangle_aspect_ratio
— Functiontriangle_aspect_ratio(inradius::Number, circumradius::Number) -> Number
Computes the aspect ratio of a triangle with inradius inradius
and circumradius circumradius
. The aspect ratio is given by
\[\tau = \dfrac{r_i}{r},\]
where $r_i$ is the inradius and $r$ is the circumradius.
triangle_aspect_ratio(p, q, r) -> Number
Computes the aspect ratio of the triangle with coordinates (p, q, r)
.
DelaunayTriangulation.triangle_radius_edge_ratio
— Functiontriangle_radius_edge_ratio(circumradius::Number, ℓmin::Number) -> Number
Computes the radius-edge ratio of a triangle with circumradius circumradius
and minimum edge length ℓmin
, given by
\[\rho = \dfrac{r}{\ell_{\min}},\]
where $r$ is the circumradius and $\ell_{\min}$ is the shortest edge length.
triangle_radius_edge_ratio(p, q, r) -> Number
Computes the radius-edge ratio of the triangle with coordinates (p, q, r)
.
DelaunayTriangulation.triangle_centroid
— Functiontriangle_centroid(p, q, r) -> (Number, Number)
Computes the centroid of a triangle with vertices p
, q
, and r
, given by
\[c = \dfrac{p + q + r}{3}.\]
DelaunayTriangulation.triangle_angles
— Functiontriangle_angles(p, q, r) -> (Number, Number, Number)
Computes the angles of a triangle with vertices p
, q
, and r
. The formula for, say, the angle at p
is given by
\[\theta_1 = \arctan\left(\dfrac{2A}{\left(p - q\right)\cdot\left(p - r\right)}\right),\]
where A
is the area of the triangle. The angles are returned in sorted order.
DelaunayTriangulation.triangle_lengths
— Functiontriangle_lengths(p, q, r) -> (Number, Number, Number)
Computes the lengths of the edges of the triangle with coordinates p
, q
, r
. The lengths are returned in sorted order.
DelaunayTriangulation.triangle_circumcenter
— Functiontriangle_circumcenter(p, q, r, A=triangle_area(p, q, r)) -> (Number, Number)
Computes the circumcenter of the triangle with coordinates (p, q, r)
. The circumcenter is given by
\[c_x = r_x + \dfrac{d_{11}d_{22} - d_{12}d_{21}}{4A}, \quad c_y = r_y + \dfrac{e_{11}e_{22} - e_{12}e_{21}}{4A},\]
where $d_{11} = \|p - r\|_2^2$, $d_{12} = p_y - r_y$, $d_{21} = \|q - r\|_2^2$, $d_{22} = q_y - r_y$, $e_{11} = p_x - r_x$ $e_{12} = d_{11}$, $e_{21} = q_x - r_x$, and $e_{22} = d_{21}$.
triangle_circumcenter(tri::Triangulation, T) -> (Number, Number)
Computes the circumcenter of the triangle T
in the triangulation tri
.
DelaunayTriangulation.triangle_offcenter
— Functiontriangle_offcenter(p, q, r, c₁=triangle_circumcenter(p, q, r), β=1.0) -> (Number, Number)
Computes the off-center of the triangle (p, q, r)
.
Arguments
p
,q
,r
: The coordinates of the triangle, given in counter-clockwise order.c₁=triangle_circumcenter(p, q, r)
: The circumcenter of the triangle.β=1.0
: The radius-edge ratio cutoff.
Output
cx
: The x-coordinate of the off-center.cy
: The y-coordinate of the off-center.
In the original paper, the off-center is defined to instead be the circumcenter if it the triangle pqc₁
has radius-edge ratio less than β
. Here, we just let the off-center be the point c
so that pqc
has radius-edge ratio of exactly β
.
DelaunayTriangulation.triangle_edge_midpoints
— Functiontriangle_edge_midpoints(p, q, r) -> (Number, Number), (Number, Number), (Number, Number)
Computes the midpoints of the edges of the triangle with coordinates (p, q, r)
.
DelaunayTriangulation.triangle_sink
— Functiontriangle_sink(tri::Triangulation, T; predicates::AbstractPredicateKernel=AdaptiveKernel()) -> (Number, Number)
Computes the sink of the triangle T
in tri
. See this paper for more information. Use the predicates
argument to control how predicates are computed.
Extended help
Sinks were introduced in this paper. For a given triangle T
, the sink of T
is defined as follows:
- If
c
, the circumcenter ofT
, is in the interior ofT
, then the sink ofT
isT
. - If
T
is a boundary triangle, then the sink ofT
isT
. - If neither 1 or 2, then the sink is defined as the sink of the triangle
V
, whereV
is the triangle adjoining the edge ofT
which intersects the linemc
, wherem
is the centroid ofT
.
In cases where the triangulation has holes, this definition can lead to loops. In such a case, we just pick one of the triangles in the loop as the sink triangle.
DelaunayTriangulation.triangle_orthocenter
— Functiontriangle_orthocenter(tri::Triangulation, T) -> NTuple{2, Number}
Finds the triangle orthocenter of T
. In particular, the point (ox, oy)
equidistant from each of the points of T
with respect to the power distance.
DelaunayTriangulation.triangle_orthoradius_squared
— Functiontriangle_orthoradius_squared(p, q, r, a, b, c) -> Number
Computes the squared orthoradius of the triangle (p, q, r)
with weights a
, b
, and c
. Note that this number may be negative.
DelaunayTriangulation.TriangulationStatistics
— TypeTriangulationStatistics{T,V,I}
A struct containing statistics about a triangulation.
Fields
num_vertices::I
: The number of vertices in the triangulation.num_solid_vertices::I
: The number of solid vertices in the triangulation.num_ghost_vertices::I
: The number of ghost vertices in the triangulation.num_edges::I
: The number of edges in the triangulation.num_solid_edges::I
: The number of solid edges in the triangulation.num_ghost_edges::I
: The number of ghost edges in the triangulation.num_triangles::I
: The number of triangles in the triangulation.num_solid_triangles::I
: The number of solid triangles in the triangulation.num_ghost_triangles::I
: The number of ghost triangles in the triangulation.num_boundary_segments::I
: The number of boundary segments in the triangulation.num_interior_segments::I
: The number of interior segments in the triangulation.num_segments::I
: The number of segments in the triangulation.num_convex_hull_vertices::I
: The number of vertices on the convex hull of the triangulation.smallest_angle::V
: The smallest angle in the triangulation.largest_angle::V
: The largest angle in the triangulation.smallest_area::V
: The smallest area of a triangle in the triangulation.largest_area::V
: The largest area of a triangle in the triangulation.smallest_radius_edge_ratio::V
: The smallest radius-edge ratio of a triangle in the triangulation.largest_radius_edge_ratio::V
: The largest radius-edge ratio of a triangle in the triangulation.area::V
: The total area of the triangulation.individual_statistics::Dict{T,IndividualTriangleStatistics{V}}
: A map from triangles in the triangulation to their individual statistics. SeeIndividualTriangleStatistics
.
Constructors
To construct these statistics, use statistics
, which you call as statistics(tri::Triangulation)
.
DelaunayTriangulation.statistics
— Functionstatistics(tri::Triangulation) -> TriangulationStatistics
Returns a TriangulationStatistics
object containing statistics about the triangulation tri
.
DelaunayTriangulation.num_boundary_segments
— Functionnum_boundary_segments(stats::TriangulationStatistics)
Returns the numboundarysegments field from the TriangulationStatistics` stats
.
DelaunayTriangulation.num_interior_segments
— Functionnum_interior_segments(stats::TriangulationStatistics)
Returns the numinteriorsegments field from the TriangulationStatistics` stats
.
DelaunayTriangulation.num_segments
— Functionnum_segments(stats::TriangulationStatistics)
Returns the num_segments field from the TriangulationStatistics` stats
.
DelaunayTriangulation.num_convex_hull_vertices
— Functionnum_convex_hull_vertices(stats::TriangulationStatistics)
Returns the numconvexhull_vertices field from the TriangulationStatistics` stats
.
DelaunayTriangulation.get_smallest_angle
— Functionget_smallest_angle(stats::TriangulationStatistics)
Returns the smallest_angle field from the TriangulationStatistics
stats
.
DelaunayTriangulation.get_largest_angle
— Functionget_largest_angle(stats::TriangulationStatistics)
Returns the largest_angle field from the TriangulationStatistics
stats
.
DelaunayTriangulation.get_smallest_area
— Functionget_smallest_area(stats::TriangulationStatistics)
Returns the smallest_area field from the TriangulationStatistics
stats
.
DelaunayTriangulation.get_largest_area
— Functionget_largest_area(stats::TriangulationStatistics)
Returns the largest_area field from the TriangulationStatistics
stats
.
DelaunayTriangulation.get_smallest_radius_edge_ratio
— Functionget_smallest_radius_edge_ratio(stats::TriangulationStatistics)
Returns the smallestradiusedge_ratio field from the TriangulationStatistics
stats
.
DelaunayTriangulation.get_largest_radius_edge_ratio
— Functionget_largest_radius_edge_ratio(stats::TriangulationStatistics)
Returns the largestradiusedge_ratio field from the TriangulationStatistics
stats
.
DelaunayTriangulation.get_individual_statistics
— Functionget_individual_statistics(stats::TriangulationStatistics)
Returns the individual_statistics field from the TriangulationStatistics
stats
.
DelaunayTriangulation.get_lengths
— Functionget_lengths(stats::TriangulationStatistics, T)
Returns the lengths field from the individual triangle statistics for the triangle T
in the TriangulationStatistics
stats
.
DelaunayTriangulation.get_circumcenter
— Functionget_circumcenter(stats::TriangulationStatistics, T)
Returns the circumcenter field from the individual triangle statistics for the triangle T
in the TriangulationStatistics
stats
.
DelaunayTriangulation.get_circumradius
— Functionget_circumradius(stats::TriangulationStatistics, T)
Returns the circumradius field from the individual triangle statistics for the triangle T
in the TriangulationStatistics
stats
.
DelaunayTriangulation.get_angles
— Functionget_angles(stats::TriangulationStatistics, T)
Returns the angles field from the individual triangle statistics for the triangle T
in the TriangulationStatistics
stats
.
DelaunayTriangulation.get_radius_edge_ratio
— Functionget_radius_edge_ratio(stats::TriangulationStatistics, T)
Returns the radiusedgeratio field from the individual triangle statistics for the triangle T
in the TriangulationStatistics
stats
.
DelaunayTriangulation.get_edge_midpoints
— Functionget_edge_midpoints(stats::TriangulationStatistics, T)
Returns the edge_midpoints field from the individual triangle statistics for the triangle T
in the TriangulationStatistics
stats
.
DelaunayTriangulation.get_aspect_ratio
— Functionget_aspect_ratio(stats::TriangulationStatistics, T)
Returns the aspect_ratio field from the individual triangle statistics for the triangle T
in the TriangulationStatistics
stats
.
DelaunayTriangulation.get_inradius
— Functionget_inradius(stats::TriangulationStatistics, T)
Returns the inradius field from the individual triangle statistics for the triangle T
in the TriangulationStatistics
stats
.
DelaunayTriangulation.get_perimeter
— Functionget_perimeter(stats::TriangulationStatistics, T)
Returns the perimeter field from the individual triangle statistics for the triangle T
in the TriangulationStatistics
stats
.
DelaunayTriangulation.get_offcenter
— Functionget_offcenter(stats::TriangulationStatistics, T)
Returns the offcenter field from the individual triangle statistics for the triangle T
in the TriangulationStatistics
stats
.
DelaunayTriangulation.get_sink
— Functionget_sink(stats::TriangulationStatistics, T)
Returns the sink field from the individual triangle statistics for the triangle T
in the TriangulationStatistics
stats
.
DelaunayTriangulation.get_minimum_angle
— Functionget_minimum_angle(stats::TriangulationStatistics, T) -> Float64
Returns the minimum angle of T
from stats
.
DelaunayTriangulation.get_maximum_angle
— Functionget_maximum_angle(stats::TriangulationStatistics, T) -> Float64
Returns the maximum angle of T
from stats
.
DelaunayTriangulation.get_median_angle
— Functionget_median_angle(stats::TriangulationStatistics, T) -> Float64
Returns the median angle of T
from stats
.
DelaunayTriangulation.get_all_stat
— Functionget_all_stat(stats::TriangulationStatistics, stat::Symbol) -> Vector
Returns a vector of the statistic stat
for each triangle in stats
.