bezier.hazmat.triangle_intersection module¶
Pure Python helper methods for bezier.triangle
.
-
bezier.hazmat.triangle_intersection.
newton_refine_solve
(jac_both, x_val, triangle_x, y_val, triangle_y)¶ Helper for
newton_refine()
.We have a system:
[A C][ds] = [E] [B D][dt] [F]
This is not a typo,
A->B->C->D
matches the Fortran-ordering of the data injac_both
. We solve directly rather than using a linear algebra utility:ds = (D E - C F) / (A D - B C) dt = (A F - B E) / (A D - B C)
- Parameters
jac_both (numpy.ndarray) – A
4 x 1
matrix of entries in a Jacobian.x_val (float) – An
x
-value we are trying to reach.triangle_x (float) – The actual
x
-value we are currently at.y_val (float) – An
y
-value we are trying to reach.triangle_y (float) – The actual
y
-value we are currently at.
- Returns
The pair of values the solve the linear system.
- Return type
-
bezier.hazmat.triangle_intersection.
newton_refine
(nodes, degree, x_val, y_val, s, t)¶ Refine a solution to \(B(s, t) = p\) using Newton’s method.
Note
There is also a Fortran implementation of this function, which will be used if it can be built.
Computes updates via
\[\begin{split}\left[\begin{array}{c} 0 \\ 0 \end{array}\right] \approx \left(B\left(s_{\ast}, t_{\ast}\right) - \left[\begin{array}{c} x \\ y \end{array}\right]\right) + \left[\begin{array}{c c} B_s\left(s_{\ast}, t_{\ast}\right) & B_t\left(s_{\ast}, t_{\ast}\right) \end{array}\right] \left[\begin{array}{c} \Delta s \\ \Delta t \end{array}\right]\end{split}\]For example, (with weights \(\lambda_1 = 1 - s - t, \lambda_2 = s, \lambda_3 = t\)) consider the triangle
\[\begin{split}B(s, t) = \left[\begin{array}{c} 0 \\ 0 \end{array}\right] \lambda_1^2 + \left[\begin{array}{c} 1 \\ 0 \end{array}\right] 2 \lambda_1 \lambda_2 + \left[\begin{array}{c} 2 \\ 0 \end{array}\right] \lambda_2^2 + \left[\begin{array}{c} 2 \\ 1 \end{array}\right] 2 \lambda_1 \lambda_3 + \left[\begin{array}{c} 2 \\ 2 \end{array}\right] 2 \lambda_2 \lambda_1 + \left[\begin{array}{c} 0 \\ 2 \end{array}\right] \lambda_3^2\end{split}\]and the point \(B\left(\frac{1}{4}, \frac{1}{2}\right) = \frac{1}{4} \left[\begin{array}{c} 5 \\ 5 \end{array}\right]\).
Starting from the wrong point \(s = \frac{1}{2}, t = \frac{1}{4}\), we have
\[\begin{split}\begin{align*} \left[\begin{array}{c} x \\ y \end{array}\right] - B\left(\frac{1}{2}, \frac{1}{4}\right) &= \frac{1}{4} \left[\begin{array}{c} -1 \\ 2 \end{array}\right] \\ DB\left(\frac{1}{2}, \frac{1}{4}\right) &= \frac{1}{2} \left[\begin{array}{c c} 3 & 2 \\ 1 & 6 \end{array}\right] \\ \Longrightarrow \left[\begin{array}{c} \Delta s \\ \Delta t \end{array}\right] &= \frac{1}{32} \left[\begin{array}{c} -10 \\ 7 \end{array}\right] \end{align*}\end{split}\]>>> nodes = np.asfortranarray([ ... [0.0, 1.0, 2.0, 2.0, 2.0, 0.0], ... [0.0, 0.0, 0.0, 1.0, 2.0, 2.0], ... ]) >>> triangle = bezier.Triangle(nodes, degree=2) >>> triangle.is_valid True >>> (x_val,), (y_val,) = triangle.evaluate_cartesian(0.25, 0.5) >>> x_val, y_val (1.25, 1.25) >>> s, t = 0.5, 0.25 >>> new_s, new_t = newton_refine(nodes, 2, x_val, y_val, s, t) >>> 32 * (new_s - s) -10.0 >>> 32 * (new_t - t) 7.0
- Parameters
nodes (numpy.ndarray) – Array of nodes in a triangle.
degree (int) – The degree of the triangle.
x_val (float) – The \(x\)-coordinate of a point on the triangle.
y_val (float) – The \(y\)-coordinate of a point on the triangle.
s (float) – Approximate \(s\)-value to be refined.
t (float) – Approximate \(t\)-value to be refined.
- Returns
The refined \(s\) and \(t\) values.
- Return type
-
bezier.hazmat.triangle_intersection.
update_locate_candidates
(candidate, next_candidates, x_val, y_val, degree)¶ Update list of candidate triangles during geometric search for a point.
Note
This is used only as a helper for
locate_point()
.Checks if the point
(x_val, y_val)
is contained in thecandidate
triangle. If not, this function does nothing. If the point is contained, the four subdivided triangles fromcandidate
are added tonext_candidates
.- Parameters
candidate (
Tuple
[float
,float
,float
,numpy.ndarray
]) –A 4-tuple describing a triangle and its centroid / width. Contains
Three times centroid
x
-valueThree times centroid
y
-value”Width” of parameter space for the triangle
Control points for the triangle
next_candidates (list) – List of “candidate” sub-triangles that may contain the point being located.
x_val (float) – The
x
-coordinate being located.y_val (float) – The
y
-coordinate being located.degree (int) – The degree of the triangle.
-
bezier.hazmat.triangle_intersection.
mean_centroid
(candidates)¶ Take the mean of all centroids in set of reference triangles.
Note
This is used only as a helper for
locate_point()
.- Parameters
candidates (
List
[Tuple
[float
,float
,float
,numpy.ndarray
] ]) –List of 4-tuples, each of which has been produced by
locate_point()
. Each 4-tuple containsThree times centroid
x
-valueThree times centroid
y
-value”Width” of a parameter space for a triangle
Control points for a triangle
We only use the first two values, which are triple the desired value so that we can put off division by three until summing in our average. We don’t use the other two values, they are just an artifact of the way
candidates
is constructed by the caller.- Returns
The mean of all centroids.
- Return type
-
bezier.hazmat.triangle_intersection.
locate_point
(nodes, degree, x_val, y_val)¶ Locate a point on a triangle.
Note
There is also a Fortran implementation of this function, which will be used if it can be built.
Does so by recursively subdividing the triangle and rejecting sub-triangles with bounding boxes that don’t contain the point. After the sub-triangles are sufficiently small, uses Newton’s method to narrow in on the pre-image of the point.
- Parameters
nodes (numpy.ndarray) – Control points for Bézier triangle (assumed to be two-dimensional).
degree (int) – The degree of the triangle.
x_val (float) – The \(x\)-coordinate of a point on the triangle.
y_val (float) – The \(y\)-coordinate of a point on the triangle.
- Returns
The \(s\) and \(t\) values corresponding to
x_val
andy_val
orNone
if the point is not on thetriangle
.- Return type
-
bezier.hazmat.triangle_intersection.
same_intersection
(intersection1, intersection2, wiggle=9.094947017729282e-13)¶ Check if two intersections are close to machine precision.
Note
This is a helper used only by
verify_duplicates()
, which in turn is only used bygeneric_intersect()
.- Parameters
intersection1 (Intersection) – The first intersection.
intersection2 (Intersection) – The second intersection.
wiggle (
Optional
[float
]) – The amount of relative error allowed in parameter values.
- Returns
Indicates if the two intersections are the same to machine precision.
- Return type
-
bezier.hazmat.triangle_intersection.
verify_duplicates
(duplicates, uniques)¶ Verify that a set of intersections had expected duplicates.
Note
This is a helper used only by
generic_intersect()
.- Parameters
duplicates (
List
[Intersection
]) – List of intersections corresponding to duplicates that were filtered out.uniques (
List
[Intersection
]) – List of “final” intersections with duplicates filtered out.
- Raises
ValueError – If the
uniques
are not actually all unique.ValueError – If one of the
duplicates
does not correspond to an intersection inuniques
.ValueError – If a duplicate occurs only once but does not have exactly one of
s
andt
equal to0.0
.ValueError – If a duplicate occurs three times but does not have exactly both
s == t == 0.0
.ValueError – If a duplicate occurs a number other than one or three times.
-
bezier.hazmat.triangle_intersection.
verify_edge_segments
(edge_infos)¶ Verify that the edge segments in an intersection are valid.
Note
This is a helper used only by
generic_intersect()
.- Parameters
edge_infos (
Optional
[list
]) – List of “edge info” lists. Each list represents a curved polygon and contains 3-tuples of edge index, start and end (see the output ofends_to_curve()
).- Raises
ValueError – If two consecutive edge segments lie on the same edge index.
ValueError – If the start and end parameter are “invalid” (they should be between 0 and 1 and start should be strictly less than end).
-
bezier.hazmat.triangle_intersection.
add_edge_end_unused
(intersection, duplicates, intersections)¶ Add intersection that is classified as “unused” and on an edge end.
This is a helper for
add_intersection()
for intersections classified asCOINCIDENT_UNUSED
. It assumes thatintersection
will have at least one ofs == 0.0
ort == 0.0
A “misclassified” intersection in
intersections
that matchesintersection
will be the “same” if it matches bothindex_first
andindex_second
and if it matches the start index exactly
- Parameters
intersection (Intersection) – An intersection to be added.
duplicates (
List
[Intersection
]) – List of duplicate intersections.intersections (
List
[Intersection
]) – List of “accepted” (i.e. non-duplicate) intersections.
-
bezier.hazmat.triangle_intersection.
check_unused
(intersection, duplicates, intersections)¶ Check if a “valid”
intersection
is already inintersections
.This assumes that
intersection
will have at least one ofs == 0.0
ort == 0.0
At least one of the intersections in
intersections
is classified asCOINCIDENT_UNUSED
.
- Parameters
intersection (Intersection) – An intersection to be added.
duplicates (
List
[Intersection
]) – List of duplicate intersections.intersections (
List
[Intersection
]) – List of “accepted” (i.e. non-duplicate) intersections.
- Returns
Indicates if the
intersection
is a duplicate.- Return type
-
bezier.hazmat.triangle_intersection.
add_intersection
(index1, s, index2, t, interior_curve, edge_nodes1, edge_nodes2, duplicates, intersections)¶ Create an
Intersection
and append.The intersection will be classified as either a duplicate or a valid intersection and appended to one of
duplicates
orintersections
depending on that classification.- Parameters
index1 (int) – The index (among 0, 1, 2) of the first edge in the intersection.
s (float) – The parameter along the first curve of the intersection.
index2 (int) – The index (among 0, 1, 2) of the second edge in the intersection.
t (float) – The parameter along the second curve of the intersection.
interior_curve (
Optional
[IntersectionClassification
]) – The classification of the intersection, if known. IfNone
, the classification will be computed.edge_nodes1 (
Tuple
[numpy.ndarray
,numpy.ndarray
,numpy.ndarray
]) – The nodes of the three edges of the first triangle being intersected.edge_nodes2 (
Tuple
[numpy.ndarray
,numpy.ndarray
,numpy.ndarray
]) – The nodes of the three edges of the second triangle being intersected.duplicates (
List
[Intersection
]) – List of duplicate intersections.intersections (
List
[Intersection
]) – List of “accepted” (i.e. non-duplicate) intersections.
-
bezier.hazmat.triangle_intersection.
classify_coincident
(st_vals, coincident)¶ Determine if coincident parameters are “unused”.
Note
This is a helper for
triangle_intersections()
.In the case that
coincident
isTrue
, then we’ll have two sets of parameters \((s_1, t_1)\) and \((s_2, t_2)\).If one of \(s_1 < s_2\) or \(t_1 < t_2\) is not satisfied, the coincident segments will be moving in opposite directions, hence don’t define an interior of an intersection.
Warning
In the “coincident” case, this assumes, but doesn’t check, that
st_vals
is2 x 2
.- Parameters
st_vals (numpy.ndarray) –
2 X N
array of intersection parameters.coincident (bool) – Flag indicating if the intersections are the endpoints of coincident segments of two curves.
- Returns
The classification of the intersections.
- Return type
-
bezier.hazmat.triangle_intersection.
should_use
(intersection)¶ Check if an intersection can be used as part of a curved polygon.
Will return
True
if the intersection is classified asFIRST
,SECOND
orCOINCIDENT
or if the intersection is classified is a corner / edge end which is classified asTANGENT_FIRST
orTANGENT_SECOND
.- Parameters
intersection (Intersection) – An intersection to be added.
- Returns
Indicating if the intersection will be used.
- Return type
-
bezier.hazmat.triangle_intersection.
triangle_intersections
(edge_nodes1, edge_nodes2, all_intersections)¶ Find all intersections among edges of two triangles.
This treats intersections which have
s == 1.0
ort == 1.0
as duplicates. The duplicates may be checked by the caller, e.g. byverify_duplicates()
.- Parameters
edge_nodes1 (
Tuple
[numpy.ndarray
,numpy.ndarray
,numpy.ndarray
]) – The nodes of the three edges of the first triangle being intersected.edge_nodes2 (
Tuple
[numpy.ndarray
,numpy.ndarray
,numpy.ndarray
]) – The nodes of the three edges of the second triangle being intersected.all_intersections (
Callable
[numpy.ndarray
numpy.ndarray
,Tuple
[numpy.ndarray
,bool
] ]) – A helper that intersects Bézier curves. Takes the nodes of each curve as input and returns an array (2 x N
) of intersections.
- Returns
4-tuple of
The actual “unique”
Intersection
-sDuplicate
Intersection
-s encountered (these will be corner intersections)Intersections that won’t be used, such as a tangent intersection along an edge
All the intersection classifications encountered
- Return type
-
bezier.hazmat.triangle_intersection.
generic_intersect
(nodes1, degree1, nodes2, degree2, verify, all_intersections)¶ Find all intersections among edges of two triangles.
This treats intersections which have
s == 1.0
ort == 1.0
as duplicates. The duplicates will be checked byverify_duplicates()
ifverify
isTrue
.- Parameters
nodes1 (numpy.ndarray) – The nodes defining the first triangle in the intersection (assumed in \(\mathbf{R}^2\)).
degree1 (int) – The degree of the triangle given by
nodes1
.nodes2 (numpy.ndarray) – The nodes defining the second triangle in the intersection (assumed in \(\mathbf{R}^2\)).
degree2 (int) – The degree of the triangle given by
nodes2
.verify (
Optional
[bool
]) – Indicates if duplicate intersections should be checked.all_intersections (
Callable
[numpy.ndarray
numpy.ndarray
,Tuple
[numpy.ndarray
,bool
] ]) – A helper that intersects Bézier curves. Takes the nodes of each curve as input and returns an array (2 x N
) of intersections.
- Returns
3-tuple of
List of “edge info” lists. Each list represents a curved polygon and contains 3-tuples of edge index, start and end (see the output of
ends_to_curve()
).”Contained” boolean. If not
None
, indicates that one of the triangles is contained in the other.The nodes of three edges of the first triangle being intersected followed by the nodes of the three edges of the second.
- Return type
-
bezier.hazmat.triangle_intersection.
geometric_intersect
(nodes1, degree1, nodes2, degree2, verify)¶ Find all intersections among edges of two triangles.
Note
There is also a Fortran implementation of this function, which will be used if it can be built.
Uses
generic_intersect()
with theGEOMETRIC
intersection strategy.- Parameters
nodes1 (numpy.ndarray) – The nodes defining the first triangle in the intersection (assumed in \(\mathbf{R}^2\)).
degree1 (int) – The degree of the triangle given by
nodes1
.nodes2 (numpy.ndarray) – The nodes defining the second triangle in the intersection (assumed in \(\mathbf{R}^2\)).
degree2 (int) – The degree of the triangle given by
nodes2
.verify (
Optional
[bool
]) – Indicates if duplicate intersections should be checked.
- Returns
3-tuple of
List of “edge info” lists. Each list represents a curved polygon and contains 3-tuples of edge index, start and end (see the output of
ends_to_curve()
).”Contained” boolean. If not
None
, indicates that one of the triangles is contained in the other.The nodes of three edges of the first triangle being intersected followed by the nodes of the three edges of the second.
- Return type
-
bezier.hazmat.triangle_intersection.
algebraic_intersect
(nodes1, degree1, nodes2, degree2, verify)¶ Find all intersections among edges of two triangles.
Uses
generic_intersect()
with theALGEBRAIC
intersection strategy.- Parameters
nodes1 (numpy.ndarray) – The nodes defining the first triangle in the intersection (assumed in \(\mathbf{R}^2\)).
degree1 (int) – The degree of the triangle given by
nodes1
.nodes2 (numpy.ndarray) – The nodes defining the second triangle in the intersection (assumed in \(\mathbf{R}^2\)).
degree2 (int) – The degree of the triangle given by
nodes2
.verify (
Optional
[bool
]) – Indicates if duplicate intersections should be checked.
- Returns
3-tuple of
List of “edge info” lists. Each list represents a curved polygon and contains 3-tuples of edge index, start and end (see the output of
ends_to_curve()
).”Contained” boolean. If not
None
, indicates that one of the triangles is contained in the other.The nodes of three edges of the first triangle being intersected followed by the nodes of the three edges of the second.
- Return type