Algorithms and tolerance settings
Some of the operations conveniently accessible via Shapes functions have non-trivial algorithms implementing them under the hood. Often, these algorithms cannot be expressed without reference to tolerance values, values that often determine tradeoffs that the user may have an interest in controlling. This part of the documentation begins by providing a guide to the various tolerace settings, and then follows descriptions of selected algorithms.
See the
man page for instructions for how to set the tolerance parameters.
Tolerance parameters
Maximum integration step length to use when integrating along paths. See also
dtmin.
Tolerance for distance computations.
Minimum integration step to use when integrating along paths. See also
arcdelta.
Tolerance used when a 3d scene is viewed in 2d.
Tolerance used when a 3d scene is viewed in 2d.
The function
..Shapes..Geometry..winding computes the winding number of a path with respect to a given point. It is only defined for closed paths. If several simple paths are grouped into a
§MultiPath, the winding number of the simple paths are added to yield the result for the combined path. The algorithm is fairly simple. Note that the result is ill-defined when the path passes through the point.
Let us first consider it for paths with only straight line segments. For this presentation, let the origin denote the point about with the winding number shall be computed, and imagine the plane separated into four quadrants meeting at the origin.
The algorithm is initialized by introducing a count which is initialized to zero, and a border between two quadrants that do not contain the starting point of the path is selected for stepping the counter.
The algorithm then traces the path, stepping from one path point to the next, keeping track of which quadrant the current point belongs to. If the quadrant does not change, the counter does not change. If the quadrant changes to an adjacent quadrant, the counter is changed appropriately if there is a transition across the selected border. If the quadrant changes to the opposite one, the step is divided into two parts by first going to the point on the segment that is closest to the origin; this will ensure that the quadrant changes to an adjacent quadrant two times, and the situation is handled accordingly.
Let us now consider paths with curved segments. The basic tool here is an algorithm to test whether the origin is included in the convex hull of the Bezier control points. If the origin is not contained, the segment may be treated as a straight line, and otherwise the segment may be sub-divided. Segments that become tiny (after many sub-dividings) indicate an ill-defined situation, and such segments are simply treated as if they were straight. The tinyness test uses a fixed tolerance that is related to the length of the curved segment before any sub-dividing.
It is important that the convex hull membership test prefers to consider the origin included, for if the origin is on the line between the first and final control point, the straight-line treatment of the segment is arbitrary, while the curved segment may be clearly on one side of the origin. That is, sub-dividing must be preferred when in doubt about curved segments. The tolerance used here is fixed, and related to the separation of the control points that generate the convex hull.
Shapes can search a point on a path which is closest to a given point or another path. The algorithms are designed to be similar in 2d and 3d to reduce development time and facilitate maintenance of the code. Hence, only when the difference is important, we shall be explicit regarding dimensions.
Since the distance between a line segment and the given point can be computed cheaply, the algorithm begins by scanning all the endpoints of the cubic spline segments and any straight line segments for an upper bound on the distance between the path and the point (recording where this distance was obtained). All curved segments are put in a job queue.
Before we continue, we just make a note of the representation used in the job queue. For each curved segment we make a base element, containing various precomputed coefficients of the spline position and its first and second derivatives. A job item is then a structure with a reference to a base element and an interval of spline times on the segment. Hence, a job item is a sub-segment of one of the original curved segments, and we shall generally speak of job item and sub-segment interchangeably..
To process the job queue, we need to compute lower bounds on the distance between sub-segment and the point, and a way to split sub-segments to obtain better estimates of the distance. With these tools at hand, a branch-and-bound algorithm is easily constructed, and we shall not go into the details of this. (Just note that each point where a sub-segment is split will be used to obtain improved upper bounds on the distance.)
To obtain lower bounds for a sub-segment, we use that the sub-segment is contained in the convex hull of the control points representation of the sub-segment. The distance between the convex hull of a set of generators and the given point can be computed in several ways, and one should use methods that are robust and efficient in the case of only four generators. Iterative improvement algorithms such as Gilbert's algorithm are avoided since it is thought that sub-optimal estimates could lead to poor convergence in ill-conditioned problems. This thought should be verified by experiments. In the meantime, we could have used the recently developed active set polytope distance algorithm (or a version of it tailored to the case of one polytope having just one generator), but we don't since it was not available at the time when the path-to-point algorithm was written. Instead, we use a scheme that enumerates the facets of the polytope, and think that this may actually be quite efficient when there is only four generators.
Computing where to split a segment is rather tricky compared to obtaining the lower bounds. The reason is that poor strategies may lead to very poor convergence (or no convergence at all) of the branch and bound algorithm. For instance, it is necessary that the splitting always result in segments that are at most some fixed factor less than one of the size of the sub-segment being split. This property will ensure that all sub-segments will eventually be shorter than
arcdelta, in which case no further splitting is made. The splitting should also occur at a local optimum along the sub-segment, as this will hopefully yield two new sub-segments which are both further away (from the given point) than the point where the split occurred. To find a good local optimizer, we first use a sparse grid of points, and from the best of these we use Newton's method.
Path-to-path approximation is a generalization of path-to-path intersections in 2d, which is both robust in the sense that no intersection is needed for the result to be well defined, and which generalizes nicely to 3d (path-to-path intersection in 3d is not a very useful concept). To generalize the intersection computation in 2d some care must be taken to insure that the first intersection is returned if there are more than one intersection, but this causes only minor thought to implement and we shall not discuss this further.
The algorithm reminds much of the path-to-point algorithm, with the major difference that the search is now performed in a two-dimensional space (one spline time along each path). The distance between pairs of straight line segments may still be computed cheaply during the initialization of the algorithm, but any other combination of straight and curved segments are pushed on the job queue.
The same type of two-layered representation of the job items is used here as in the path-to-point algorithm, with one base element for each pair of original spline segments, and the two major tasks are the computation of lower bounds and splitting times. As in the path-to-point algorithm, lower bounds are obtained by computing the distance between the convex hulls of the control point representation of the two sub-segments in a job item. For this, an active set algorithm is used. The active set algorithm deserves its own documentation, but it remains to be written (at this time, we just remark that
disttol is related to Lagrange multipliers by a nice geometric interpretation of a scaled gradient). When computing where to split a job item, it is important that both sub-segments are divided so that they both eventually become smaller than
arcdelta. (Just making the “area” spanned by the job item's two dimensions tend to zero is not enough, as there is no efficient way to estimate the distance between a long and a very short sub-segment — two very short segments, on the other hand, are well approximated by their control point convex hull distance.) Again (as in the path-to-point algorithm) we first use a grid and then Newton's method.