Guarantees Concerning Geometric Objects with Imprecise Points by Jeff Sember B.Sc., School of Computing Science, Simon Fraser University, 2004 M.Sc., School of Computing Science, Simon Fraser University, 2006 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF Doctor of Philosophy in THE FACULTY OF GRADUATE STUDIES (Computer Science) The University Of British Columbia (Vancouver) October 2011 c Jeff Sember, 2011 Abstract Traditional geometric algorithms are often presented as if input imprecision does not exist, even though it is often unavoidable. We argue that in some cases, it may be desirable for geometric algorithms to treat this imprecision as an explicit component of the input, and to reflect this imprecision in the output. Starting with three problems from computational geometry whose inputs are planar point sets (Voronoi diagrams, convex hulls, and smallest bounding discs), we recast these as problems where each input point’s location is imprecise, but known to lie within a particular region of uncertainty. Where algorithms to solve each of the original problems produce a single geometric object as output, the algorithms that we present typically produce either guaranteed or possible output objects. A guaranteed object represents qualities that can be guaranteed for every point set that is consistent with the uncertain regions, and a possible object represents qualities that exist for at least one such point set. By dealing with input imprecision explicitly, these guaranteed and possible objects can represent a more accurate view of what can be reliably inferred from the input data. ii Preface Portions of Chapter 2 have been published: Sember, J. and Evans, W. Guaranteed Voronoi diagrams of uncertain sites. Proceedings of the 20th Canadian Conference on Computational Geometry, pages 207–210, 2008. The research and writing is mine, except for the proof of Lemma 3, which was contributed by Will Evans. Will also assisted in proofreading the paper. iii Table of Contents Abstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ii Preface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iii Table of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iv List of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vii List of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viii Acknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xii Dedication . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiii 1 2 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.1 Uncertainty, Precision, and Accuracy . . . . . . . . . . . . . . . . 1 1.2 Robustness Issues . . . . . . . . . . . . . . . . . . . . . . . . . . 3 1.3 Geometric and Topological Concepts . . . . . . . . . . . . . . . . 5 1.4 Regions of Uncertainty . . . . . . . . . . . . . . . . . . . . . . . 7 1.5 Contributions of this Thesis . . . . . . . . . . . . . . . . . . . . . 9 Voronoi Diagram of Imprecise Points . . . . . . . . . . . . . . . . . 11 2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 2.1.1 Related Work . . . . . . . . . . . . . . . . . . . . . . . . 14 2.1.2 Contributions . . . . . . . . . . . . . . . . . . . . . . . . 17 Guaranteed Voronoi Diagram . . . . . . . . . . . . . . . . . . . . 17 2.2.1 17 2.2 Properties . . . . . . . . . . . . . . . . . . . . . . . . . . iv 2.2.2 A Mapping between Boundaries of Voronoi Cells . . . . . 20 2.2.3 Uncertain Discs . . . . . . . . . . . . . . . . . . . . . . . 22 2.2.4 Uncertain Polygons . . . . . . . . . . . . . . . . . . . . . 26 2.3 Possible Voronoi Diagram . . . . . . . . . . . . . . . . . . . . . 37 2.4 Guaranteed Delaunay Edges . . . . . . . . . . . . . . . . . . . . 44 2.5 Guaranteed Gabriel Edges of Uncertain Discs . . . . . . . . . . . 51 2.5.1 Gabriel Discs and Zones . . . . . . . . . . . . . . . . . . 51 2.5.2 Finding Guaranteed Gabriel Edges . . . . . . . . . . . . . 59 Closing Remarks . . . . . . . . . . . . . . . . . . . . . . . . . . 66 Convex Hull of Imprecise Points . . . . . . . . . . . . . . . . . . . . 70 3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70 3.1.1 Related Work . . . . . . . . . . . . . . . . . . . . . . . . 73 3.1.2 Contributions . . . . . . . . . . . . . . . . . . . . . . . . 74 3.2 Guaranteed Hull . . . . . . . . . . . . . . . . . . . . . . . . . . . 75 3.3 Guaranteed Hull Complexity . . . . . . . . . . . . . . . . . . . . 81 3.4 Guaranteed Hull of Aligned Similar Uncertain Regions . . . . . . 86 3.4.1 Algorithm 1 . . . . . . . . . . . . . . . . . . . . . . . . . 86 3.4.2 Algorithm 2 . . . . . . . . . . . . . . . . . . . . . . . . . 93 2.6 3 3.4.3 3.5 3.6 4 Guaranteed Hulls in R3 . . . . . . . . . . . . . . . . . . . Possible Hull . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100 3.5.1 Point and Polygon . . . . . . . . . . . . . . . . . . . . . 108 3.5.2 Two Polygons . . . . . . . . . . . . . . . . . . . . . . . . 115 3.5.3 Multiple Polygons . . . . . . . . . . . . . . . . . . . . . 124 Closing Remarks . . . . . . . . . . . . . . . . . . . . . . . . . . 124 Smallest Bounding Disc of Imprecise Points 4.1 97 . . . . . . . . . . . . . 127 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 127 4.1.1 Related Work . . . . . . . . . . . . . . . . . . . . . . . . 129 4.1.2 Contributions . . . . . . . . . . . . . . . . . . . . . . . . 130 4.2 Possible 1-Centers . . . . . . . . . . . . . . . . . . . . . . . . . 130 4.3 Query Algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . 139 4.3.1 Uncertain Axis-Aligned Rectangles . . . . . . . . . . . . 150 v 4.4 Construction Algorithm . . . . . . . . . . . . . . . . . . . . . . . 152 4.4.1 4.5 5 Uncertain Axis-Aligned Rectangles . . . . . . . . . . . . 154 Closing Remarks . . . . . . . . . . . . . . . . . . . . . . . . . . 159 Conclusion and Future Research . . . . . . . . . . . . . . . . . . . . 161 Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 163 Appendix A Proof of Lemma 4.16 . . . . . . . . . . . . . . . . . . . . 170 vi List of Tables Table A.1 Proof of Lemma 4.16, Subcases. . . . . . . . . . . . . . . . . 172 Table A.2 Proof of Lemma 4.16, Cases 1 . . . 5. . . . . . . . . . . . . . . . 172 Table A.3 Proof of Lemma 4.16, Cases 6 . . . 13. . . . . . . . . . . . . . . 173 Table A.4 Proof of Lemma 4.16, Cases 14, 15. . . . . . . . . . . . . . . . 174 vii List of Figures Figure 1.1 Valid pairs Di , Dj ∈ D. . . . . . . . . . . . . . . . . . . . . 8 Figure 1.2 Invalid pairs Di , Dj ∈ D. . . . . . . . . . . . . . . . . . . . 8 Figure 1.3 Valid (Di ) and invalid (Dj , Dk ) uncertain polygons. . . . . . 8 Figure 2.1 Voronoi diagram of a set of points. . . . . . . . . . . . . . . . 12 Figure 2.2 Guaranteed Voronoi diagram for uncertain discs. . . . . . . . 13 Figure 2.3 Delaunay triangulation of a set of points, the dual of the Voro13 Figure 2.4 noi diagram (circumcircle for one triangle is shown). . . . . . ei . . . . . . . . . . . . . . . . . . . . . . . Lemma 2.1: p ∈ R Figure 2.5 Lemma 2.4. . . . . . . . . . . . . . . . . . . . . . . . . . . . 20 Figure 2.6 Definition of δ(p). . . . . . . . . . . . . . . . . . . . . . . . 21 Figure 2.7 Lemma 2.6. . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 Figure 2.8 Lemma 2.6: pδ(p) and qδ(q) cannot intersect. . . . . . . . . . ei . . . . . . . . . Multiple edges of hi, ji on the boundary of R 23 24 Figure 2.10 Guaranteed Voronoi diagram of uncertain polygons. . . . . . 27 Figure 2.11 Lemma 2.11. . . . . . . . . . . . . . . . . . . . . . . . . . . ei . Also shown are farthest point Figure 2.12 Finding point on boundary of R 30 Figure 2.9 19 Voronoi diagram of Vi (solid lines) and bisector of u and w (dashed line). . . . . . . . . . . . . . . . . . . . . . . . . . . 31 Figure 2.13 Lemma 2.13. . . . . . . . . . . . . . . . . . . . . . . . . . . 32 Figure 2.14 Events that can occur while growing Cs . . . . . . . . . . . . . 33 ei (Cp0 equals hv, v 0 , ei in this example). 35 Figure 2.15 v 0 cannot contribute to R Figure 2.16 Possible Voronoi diagram. . . . . . . . . . . . . . . . . . . . 38 Figure 2.17 The possible cell for a site. . . . . . . . . . . . . . . . . . . . 38 Figure 2.18 Possible cell vertices lie on edges of V viii m (D). . . . . . . . . . 40 Figure 2.19 Ve {} (D) complexity for discs is Ω(n3 ) (some edges omitted for clarity). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43 Figure 2.20 Cocircular points with two valid Delaunay triangulations (Voronoi diagram drawn with dotted lines). . . . . . . . . . . . . . 44 Figure 2.21 Lemma 2.22. . . . . . . . . . . . . . . . . . . . . . . . . . . 45 Figure 2.22 Lemma 2.23. . . . . . . . . . . . . . . . . . . . . . . . . . . 46 Figure 2.23 Lemma 2.23 does not apply for nonconvex regions. . . . . . . 47 Figure 2.24 The Gabriel zone for a pair of uncertain discs. . . . . . . . . . 52 Figure 2.25 Lemma 2.31. . . . . . . . . . . . . . . . . . . . . . . . . . . 53 Figure 2.26 Lemma 2.32. . . . . . . . . . . . . . . . . . . . . . . . . . . 54 Figure 2.27 Lemma 2.33: there are exactly two transition points. . . . . . 55 Figure 2.28 Lemma 2.33: boundary points above the x-axis are left-oriented. 56 Figure 2.29 Deriving an expression for boundary points of Za,b . . . . . . . 56 Figure 2.30 Determining if Dk intersects Za,b . . . . . . . . . . . . . . . . 58 Figure 2.31 Dk intersects Gabriel disc Cp , yet is not a member of N . . . . 61 Figure 2.32 Lemma 2.35. . . . . . . . . . . . . . . . . . . . . . . . . . . 62 Figure 2.33 p1 is on the boundary of R{a,b} . . . . . . . . . . . . . . . . . 63 Figure 2.34 p1 lies on hha, bii. . . . . . . . . . . . . . . . . . . . . . . . . 64 Figure 2.35 Guaranteed spanning edges. . . . . . . . . . . . . . . . . . . 68 Figure 2.36 Applet for Chapter 2. . . . . . . . . . . . . . . . . . . . . . . 69 Figure 3.1 The convex hull of a set of points . . . . . . . . . . . . . . . 71 Figure 3.2 Convex hull, guaranteed hull, and possible hull of uncertain regions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73 Figure 3.3 L is supported by Di , and invalidated by Dj . . . . . . . . . . 76 Figure 3.4 Bitangents . . . . . . . . . . . . . . . . . . . . . . . . . . . . 77 Figure 3.5 Theorem 3.1. . . . . . . . . . . . . . . . . . . . . . . . . . . 79 Figure 3.6 Lemma 3.2. . . . . . . . . . . . . . . . . . . . . . . . . . . . 79 Figure 3.7 Lemma 3.5 . . . . . . . . . . . . . . . . . . . . . . . . . . . 82 Figure 3.8 Theorem 3.6 . . . . . . . . . . . . . . . . . . . . . . . . . . 83 Figure 3.9 Lemma 3.7: 2n − 2 hull bitangents . . . . . . . . . . . . . . 84 Figure 3.10 Some hull bitangents of a set of uncertain triangles . . . . . . 85 Figure 3.11 Hull bitangents of discs . . . . . . . . . . . . . . . . . . . . . 87 ix Figure 3.12 Bxy cannot be the stack’s topmost element. . . . . . . . . . . 97 Figure 3.13 Lemma 3.20 . . . . . . . . . . . . . . . . . . . . . . . . . . 103 Figure 3.14 Theorem 3.22 . . . . . . . . . . . . . . . . . . . . . . . . . . 104 Figure 3.15 Possible hull of uncertain polygons . . . . . . . . . . . . . . 105 Figure 3.16 Theorem 3.24 . . . . . . . . . . . . . . . . . . . . . . . . . . 107 Figure 3.17 Uncertain polygons with n vertices generating P H with n vertices . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108 Figure 3.18 Possible hull of polygon P and point s . . . . . . . . . . . . . 108 Figure 3.19 Expansion step, case (i), before and after changes to hull. . . . 110 Figure 3.20 Expansion step, case (ii), before and after changes to H. . . . 110 Figure 3.21 Expansion step, case (iii), before and after changes to H. . . . 111 Figure 3.22 Interior step, case (ii). . . . . . . . . . . . . . . . . . . . . . 111 Figure 3.23 Lemma 3.25. . . . . . . . . . . . . . . . . . . . . . . . . . . 114 Figure 3.24 The interior of uv cannot contain two new vertices of H. . . . 115 Figure 3.25 Initial polygon H; pocket lids shown as dashed lines. . . . . . 116 Figure 3.26 Hull contraction: pocket lid ai aj has been replaced by (ai . . . aj ) of J. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117 Figure 3.27 Hull expansion (ccw direction). . . . . . . . . . . . . . . . . 118 Figure 3.28 Lemma 3.28. . . . . . . . . . . . . . . . . . . . . . . . . . . 120 Figure 3.29 Expanding H remains within P H. . . . . . . . . . . . . . . . 121 Figure 3.30 Wedges induced by pocket (ai . . . aj ). . . . . . . . . . . . . . 122 Figure 3.31 Lemma 3.28 will not hold. . . . . . . . . . . . . . . . . . . . 123 Figure 3.32 Applet for Chapter 3. . . . . . . . . . . . . . . . . . . . . . . 125 Figure 4.1 Smallest bounding discs of point sets. . . . . . . . . . . . . . 129 Figure 4.2 Lemma 4.2. . . . . . . . . . . . . . . . . . . . . . . . . . . . 132 Figure 4.3 R is locally inside-tangent to C at p. . . . . . . . . . . . . . . 132 Figure 4.4 Point p is on the type III curve of Df , Di , and Dj . . . . . . . 133 Figure 4.5 Polar angle ranges associated with qi ∈ Di . . . . . . . . . . . 134 Figure 4.6 Right(i) ∩ Right(j) 6= ∅. . . . . . . . . . . . . . . . . . . . . 135 Figure 4.7 In(i) ∩ Out(j) 6= ∅. . . . . . . . . . . . . . . . . . . . . . . . 136 Figure 4.8 Possible 1-centers of uncertain discs. . . . . . . . . . . . . . 137 Figure 4.9 Parameterization of IIIf ib . . . . . . . . . . . . . . . . . . . . 139 x Figure 4.10 Terminology related to query point q. . . . . . . . . . . . . . 140 Figure 4.11 Examples of gaps at disc boundary points. . . . . . . . . . . . 141 Figure 4.12 Lemma 4.11 does not hold for u, u e. . . . . . . . . . . . . . . 143 Figure 4.13 Disc curves of left half gap function for pr . . . . . . . . . . . 144 Figure 4.14 Lemma 4.12. . . . . . . . . . . . . . . . . . . . . . . . . . . 146 Figure 4.15 Lemma 4.14. . . . . . . . . . . . . . . . . . . . . . . . . . . 148 Figure 4.16 Boundary edges for axis-aligned rectangles. . . . . . . . . . . 150 Figure 4.17 Lemma 4.18. . . . . . . . . . . . . . . . . . . . . . . . . . . 152 Figure 4.18 Possible 1-centers of uncertain rectangles. . . . . . . . . . . . 155 Figure 4.19 Lemma 4.22. . . . . . . . . . . . . . . . . . . . . . . . . . . 156 Figure 4.20 IIIf ij is a line. . . . . . . . . . . . . . . . . . . . . . . . . . 157 Figure 4.21 IIIf ij is a parabola. . . . . . . . . . . . . . . . . . . . . . . 158 Figure 4.22 Possible bounding disc, guaranteed bounding disc (uncertain discs). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 159 Figure 4.23 Possible bounding disc, guaranteed bounding disc (uncertain axis-aligned rectangles). . . . . . . . . . . . . . . . . . . . . 159 Figure 4.24 Applet for Chapter 4. . . . . . . . . . . . . . . . . . . . . . . 160 Figure A.1 Fact (iv). . . . . . . . . . . . . . . . . . . . . . . . . . . . . 171 xi Acknowledgments I am very grateful to the members of my supervisory committee. I consider myself very fortunate to have had Will Evans as my senior supervisor; his knowledge and guidance has been invaluable during my time at the University of British Columbia. I am also indebted to David Kirkpatrick, whose help during the final stages of my thesis writing has significantly improved my dissertation. Some other members of the computational geometry community that I would be remiss in not giving thanks to include Maarten Löffler, whom I had the pleasure of having a number of fruitful discussions with during his visits to Vancouver; and Binay Bhattacharya, my supervisor at Simon Fraser University, whose ongoing friendship and encouragement I value. I am grateful for the recognition and considerable financial support that I have received from the Natural Sciences and Engineering Research Council of Canada. The funding they have provided me in my graduate studies has been a key motivation for me continuing to pursue my academic goals, and I feel fortunate to live in a country that has such an institution. My PhD studies have been made more enjoyable by my friends and colleagues in the Beta lab, as well as other computer science students at UBC. Inevitably, some of these people have moved on to bigger and better things, but I hope to keep them as friends for a long time to come. Finally, I would like to thank my friends and family, for their love, encouragment, and cooking. xii Dedication For dad xiii Chapter 1 Introduction Imprecision is usually unavoidable in geometric algorithms. Despite this fact, such algorithms are typically presented as if imprecision does not exist: they implicitly assume that the input values exactly represent some true value, when they actually include some degree of imprecision. This thesis argues that geometric algorithms should treat this imprecision as an explicit component of the input, and reflect this imprecision in the output. We address some fundamental problems from computational geometry involving (planar) point inputs, and recast them into problems where each point’s precise location is unknown but is guaranteed to lie within a known region of uncertainty. Typically, these problems will come in one of two flavours: a guaranteed problem, whose solution expresses qualities that are shared by the solutions for every possible choice of precise inputs; and a possible problem, whose solution expresses qualities that are shared by the solution for at least one possible choice of precise inputs. By structuring the problems in this way, we ensure that the solutions reflect the strongest possible guarantees that can be made given the imprecision in the input. 1.1 Uncertainty, Precision, and Accuracy Geometric algorithms rarely produce answers that can be viewed with complete certainty. Three sources of imprecision that contribute to this uncertainty have been identified by Löffler [45]. In this section, we examine each of these three 1 types of imprecision, and discuss how they will be addressed in this thesis. We can label the first of the three types as measurement imprecision. When the input to a geometric algorithm is derived from objects in the physical world, there is inevitably some device or instrument (which may incorporate one or more human senses) that is used to take a measurement. Since these are real-world objects, they necessarily produce readings that are accompanied by some degree of imprecision. For example, the distance between two entrances to a termite hill may be measured with a ruler that is calibrated in sixteenths of an inch, and the location of the termite hill may be read from a GPS (Global Positioning Satellite) device that reports locations using some number of significant digits. Both instruments incorporate some degree of imprecision. The second type of imprecision, modeling imprecision, concerns the way a real-world object is represented within a geometric algorithm. The rim of a canyon cannot be described exactly with limited resources, since this would entail measuring the location of the infinitely many points defining its boundary. Instead, some finite number of samples are usually taken, and the canyon rim is then represented as a geometric curve. Such a curve is clearly only an approximation of the true canyon rim. As another example, consider the task of recording the locations of every tree on a particular farmer’s field. It is natural to construct a two-dimensional map of the trees. Such a map can only be an approximation to what is actually a portion of a three-dimensional object, the surface of the Earth. While it is true that the error introduced in this way is probably minimal (unless the field is the size of Belgium), the fact remains that in this example, some modeling imprecision cannot be avoided. The third type of imprecision concerns the way computers store and manipulate numbers to calculate results, and can be labelled computational imprecision. It is usually implicitly assumed that the inputs to an algorithm are real numbers, yet such numbers cannot be exactly represented by computers (except for some subsets of real numbers, such as integers or rational numbers). Geometric algorithms are often presented as if they are to be run on a theoretical computer that can perform exact calculations on real numbers (one such machine is the Real RAM, or Random Access Machine). In practice, real numbers must be represented and manipulated using finite resources. This is usually done by approximating a real number with a 2 floating point value, which uses (typically 32 or 64) bits to represent a finite subset of the real numbers. This can result in programs that are nonrobust: for some valid input, they produce an incorrect or inaccurate result, or fail to return a result at all (e.g., by entering an infinite loop, or crashing). An aspect of uncertainty that is closely related to precision is accuracy. Whereas the precision associated with a value determines how the value is stored (e.g., the number of digits to the right of the decimal point, or the resolution of a bitmap), accuracy represents how close the stored value is to the true value. Accuracy and precision are independent notions. For example, the number 1.9998 can be accurately, yet imprecisely, represented by the integer 2, and it can be inaccurately, yet precisely, represented by the number 3.1416. Clearly, precision without accuracy is not very useful; hence, when possible, the precision of a value (e.g. a reading, a measurement, or a calculation) is often reduced so that the value’s inaccuracy matches its imprecision. 1.2 Robustness Issues It is easy to find examples of geometric algorithms that lead to nonrobust programs. Consider the following problem: given nonparallel lines L1 and L2 in the plane, first calculate their point of intersection p = L1 ∩ L2 , then determine whether p lies on L1 . Clearly, for any nonparallel input lines, the algorithm should return true; yet when implemented using standard floating point libraries, the resulting program returns false nearly half the time [34]. Many techniques have been developed to address the robustness issues that arise when implementing geometric algorithms. Perhaps the most common method is ‘epsilon tweaking’: when deciding whether the result x of a computation is zero, the program returns true iff x is within some small constant of zero (i.e, iff |x| < ). This is usually implemented in an ad-hoc fashion, often with many different values of epsilon spread throughout a program. While this technique has the advantage of being simple, it does not solve the problem: the probability of false negatives is reduced only at the expense of increased probability of false positives, and the possibility of endless loops and program crashes still exists. A more systematic approach for addressing robustness issues is the exact com- 3 putation paradigm [77], in which a program represents real numbers exactly. This technique is not without its drawbacks: first, the numbers that are represented are restricted to a small subset of real numbers (e.g., rational numbers); and second, the use of such methods often results in much longer execution times. Another technique is epsilon geometry [30], in which standard geometric objects (e.g., points) are effectively replaced by ‘fat’ geometric objects (e.g., discs). Epsilon geometry employs interval arithmetic [52] to replace yes/no predicates with partial yes/no/maybe predicates. While robustness issues often result from rounding errors and other shortcomings inherent in standard floating point libraries, they can also be a consequence of degeneracies in the inputs. For example, geometric algorithms often assume that the intersection between a pair of line segments is either the empty set, or a single point interior to both segments. When given a pair of segments that are collinear, or whose intersection consists of one of the segment’s endpoints, the algorithm can misbehave in unforseen ways. For this reason, many geometric algorithms (including most of the ones presented in this thesis) make various general position assumptions: no two lines are collinear, no three lines intersect at a common point, and so on. Modifying an algorithm to handle such degeneracies can be complicated, and often requires significant work. A common method for dealing with degenerate inputs is to modify the inputs very slightly and at random, so that with high probability the new inputs are degeneracy-free. An implicit assumption of such perturbation schemes is that the output for the modified inputs is sufficiently close, in some sense, to that of the original inputs. In some situations, it may be possible to mitigate computational imprecision at the expense of an increase in measurement imprecision. Consider an algorithm that manipulates axis-aligned rectangles. Perhaps the calculations involved will be much simpler, or at least more robust, if the coordinates of the sides of the rectangles are integers. If this is the case, we can expand each input rectangle so that it contains the original rectangle, and has integer coordinates. If the resolution of the integer grid is fine enough, the increase in measurement imprecision will be a small price to pay for the reduction in computational imprecision. 4 1.3 Geometric and Topological Concepts In this section, we introduce some terminology concerning geometric and topological concepts that we will be using throughout the thesis. All geometric objects that we will discuss are assumed to be two-dimensional: that is, subsets of R × R. All distances will be assumed to be Euclidean, and we will denote the distance between points a and b as d(a, b). We may use the wildcard character ‘∗’ to denote a location within an expression where any appropriate index or variable can be substituted. If two or more wildcards appear in one expression, different values can be substituted for each occurrence. The polar angle of a directed line is the angle that the line makes with the positive x-axis. The polar angle of a ray (or directed segment) is the polar angle of the directed line that contains, and is in the same direction as, the ray (or segment). When we refer to the left or right side of a particular segment that is denoted by ab, we will assume the segment is directed from a to b. The polar angle of a directed line L is denoted θ(L). We will make a distinction between circles and discs, by considering a circle to be the boundary of a disc. In other words, a disc includes the points it surrounds, while a circle does not. In general, a region need not contain its boundary. If it does, then it is closed. For example, the set of points {(x, y) | 0 ≤ x ≤ 2, 1 ≤ y ≤ 3}, which represents a rectangle, is a closed region. By contrast, the set of points {(x, y) | x2 + y 2 < 9}, which represents the interior of a disc, is not closed. The closure of a region is the union of a region and its boundary. A region is bounded if it is a subset of some disc of finite radius. In other words, there exists some constant k such that no two points of the region are separated by a distance greater than k. A region that is both closed and bounded is compact. Let A and B be two regions. A and B are disjoint if they do not intersect, i.e., if A ∩ B = ∅. If A contains B, and no part of B intersects the boundary of A, then we say that B is interior to, or lies in the interior of, A. If neither of two regions is interior to the other, then we say that the regions are partially disjoint (regions that are disjoint are also considered to be partially disjoint). We will say that B 5 penetrates A if B intersects the interior of A. To summarize, in Figure 1.1, the regions on the left are both disjoint and partially disjoint; the regions in the center are only partially disjoint; and to the right, Dj is interior to Di . Regions Di and Dj penetrate each other in all but the leftmost sets of both Figure 1.1 and Figure 1.2. B is inside-tangent to A (or A has inside-tangent B) if B ⊆ A and the boundary of A intersects B. We will say that A is tangent to B at the points where the boundaries intersect. B is outside-tangent to A (or A has outside-tangent B) if B ∩ A is a non-empty subset of the boundary of A. Again, A is tangent to B at the points in this intersection. Note that in contrast to inside-tangency, outsidetangency is reflexive: B is outside-tangent to A iff A is outside-tangent to B. We will say that directed line L is right-tangent to A if A lies in the right half-plane of L, and A intersects L (we may also say that A is right-tangent to L). We define left-tangency analogously. We will use cw and ccw as abbreviations for clockwise and counterclockwise respectively. A ccw range [θ1 θ2 ] is the set of angles from θ1 ccw through θ2 . We will denote open or half-open sets as (θ1 θ2 ), [θ1 θ2 ), etc. A polygonal chain is a sequence of directed line segments {p1 p2 , p2 p3 , . . . , pn−1 pn }. For succinctness, we usually describe such a chain as a sequence of vertices {p1 , . . . , pn }. A polygonal chain is simple if it has no self-intersections.A polygon is a region whose boundary is a polygonal chain, and it is simple if its boundary is simple. We will assume that simple polygons have a ccw winding: the interior of the polygon lies to the left of each of its boundary’s directed segments. A nonsimple polygon can be problematic, since it may not have a winding direction or interior that satisfies this definition. See, for example, polygon Dj of Figure 1.3. Davenport-Schinzel sequences arise often in complexity proofs of geometric structures ([17],[67]). These are sequences of elements in which the number of times a particular pair of elements can appear in alternating sequence is limited. Formally, a sequence S is a Davenport-Schinzel sequence of order k if the following conditions hold: (i) there are at most n distinct elements in S; (ii) no two adjacent elements of S are the same; and (iii) S contains no subsequence (a, b, a, b, . . .) of length k + 2. 6 It can be shown that the maximum size of a Davenport-Schinzel sequence of order 2 is 2n − 1. For k = 3, the bound increases to nα(n), where α(·) is the extremely slowly-growing inverse Ackermann function (for all practical values of n, α(n) < 5). If, in addition to satisfying the three conditions above, S has the property that its first and last elements are distinct (or S has fewer than two elements), then S is a Davenport-Schinzel cyclic sequence of order k. 1.4 Regions of Uncertainty A region of uncertainty is a subset of the plane that is associated with a point site whose exact location is unknown, but is guaranteed to lie within the region of uncertainty. We may refer to such a region as an uncertain region, which is perhaps misleading, since the region itself is known precisely; rather, it is the location of the point within the region that is uncertain. Discs and squares are perhaps the two most common shapes for uncertain regions to take. As observed in [47], an uncertain disc is likely to arise when the coordinates of a point have been measured from some physical device (e.g., a scanner, or a GPS reading), while an uncertain square usually results when the point’s coordinates have been stored as a pair of floating point values, with each value having independent amounts of imprecision. In addition to discs and squares, in this thesis we will also be considering other shapes, such as axis-aligned rectangles and polygons. We will denote a set of uncertain regions {D1 , . . . , Dn } by D. We will assume that each Di ∈ D is compact. In addition, we will require the uncertain regions to satisfy (i) the general position assumption that no line contains more than two vertices of D, and (ii) a valid intersection property: intersections between the boundaries of distinct regions of D must be proper intersections. This second condition can be stated formally as follows: if p is a point where the boundaries of distinct regions Di and Dj of D intersect, then p must lie on the boundary separating two nonzero-area subsets of Di , where one subset is strictly inside Dj , and the other is strictly outside Dj . Figure 1.1 displays examples of pairs of regions that satisfy this condition, while Figure 1.2 displays examples of pairs of regions that do not. 7 Dj Di Di Dj Di Dj Figure 1.1: Valid pairs Di , Dj ∈ D. Di Di Di Dj Dj Dj Figure 1.2: Invalid pairs Di , Dj ∈ D. We will refer to a set of precise, pairwise distinct points drawn from uncertain regions as a feasible set. If an uncertain region Di is a disc, we will denote its origin by oi , and its radius by ri . If an uncertain region Di is a polygon, we will assume that it is simple, has a ccw winding, and has at least three vertices. Figure 1.3 illustrates some valid and invalid uncertain polygons. Di Dj Dk Figure 1.3: Valid (Di ) and invalid (Dj , Dk ) uncertain polygons. 8 1.5 Contributions of this Thesis This thesis is concerned with problems involving sets of points whose locations are uncertain, but known to lie within particular regions. We will present algorithms that when given such regions as input, produce either guaranteed or possible output objects. Together, these two objects allow us to make the strongest possible guarantees about arbitrary points in the plane, with respect to the input points. For example, consider the problem of determining the smallest axis-aligned rectangle that contains a set of input points (i.e., their ‘smallest bounding rectangle’). We can solve this simple problem in linear time, and we can be confident that for any set of input points, the rectangle produced is the correct answer. Now consider the situation when the location of each input point is imprecise, and may lie anywhere within a particular disc. We cannot hope to produce the actual smallest bounding rectangle of the points, since we do not know exactly where the points lie. The strongest guarantee we can make about any particular point in the plane is to assign it to one of three sets: points that lie within the smallest bounding rectangle of every feasible set; points that lie within the smallest bounding rectangle of at least one feasible set; and a third set, which contains any points not in either of the previous two. The guaranteed object produced by our algorithm would be the first of these sets, and the possible object would be the union of the first and second sets. Throughout this thesis, we will assume that the uncertain regions themselves are known exactly, and can be represented precisely (i.e., as a disc with known origin and radius, or as a polygon with known vertices). Since these regions are clearly mathematical entities, we will avoid introducing any modeling imprecision. We will also largely avoid addressing issues of computational imprecision, and will usually assume the use of the aforementioned Real RAM as the theoretical computing device. In addition, we will ignore any issues concerning the accuracy of the values manipulated by the algorithms. In each of the next three chapters, we address some fundamental problems from computational geometry involving planar point sites, where each site lies within a region of uncertainty. In Chapter 2, we focus on Voronoi diagrams, and two related structures: Delaunay triangulations, and Gabriel graphs. We provide complexity 9 bounds and algorithms to construct guaranteed Voronoi diagrams, and to identify guaranteed Delaunay and Gabriel edges. We introduce a new geometric structure, the possible Voronoi diagram, and give a complexity result and algorithm to construct this figure. In Chapter 3, we consider the problem of convex hulls of imprecise points. We give new results on the complexity of guaranteed hulls, and new algorithms for their construction. We also investigate the possible hull of imprecise points. This is a novel problem for uncertain regions that are nonconvex (otherwise, the possible hull is just the convex hull of the regions), and we present an algorithm that constructs the possible hull of nonconvex polygons. In Chapter 4, we address the smallest bounding discs of point sets. We describe a new geometric figure, the possible 1-centers of uncertain regions, and for uncertain discs and uncertain axis-aligned rectangles, we provide two algorithms related to this figure. We also highlight a connection between this figure and guaranteed hulls, one that can be exploited to optimize algorithms involving possible 1-centers. We conclude with a brief summary of our results, and some directions for future research. At the end of each chapter, we will provide a link to a website containing a program that illustrates some of the material presented in the chapter. Each of these programs is a Java applet which enables users to create and manipulate files of geometric objects (e.g., discs and polygons), and demonstrates some of the algorithms presented in this thesis. 10 Chapter 2 Voronoi Diagram of Imprecise Points 2.1 Introduction Voronoi diagrams are a fundamental data structure in computational geometry. Given a set of points {p1 , . . . , pn }, the Voronoi diagram (Figure 2.1) partitions the plane into n regions, or cells, such that each point in cell i is closer to pi than to any p6=i . One application of Voronoi diagrams is the post-office problem [42], in which the sites correspond to post office locations, and the Voronoi diagram indicates which post office is closest to a particular query point. Additional data structures can be used in conjunction with Voronoi diagrams to answer such queries efficiently (typically in logarithmic time). There are many algorithms for constructing Voronoi diagrams. Techniques include the divide and conquer approach employed by Shamos and Hoey [65], Fortune’s plane sweep method [26], and the higher dimensional embedding of Brown [12], which shows that the Voronoi diagram of a set of (planar) points can be recovered from a 3-dimensional hull induced by the points. Lee and Schachter [43] generate Voronoi diagrams by first constructing their dual, the Delaunay triangulation. Each of these algorithms runs in O(n log n) time, and (since Voronoi di- 11 agrams can be used to sort a set of numbers) is therefore worst-case optimal.1 Voronoi diagrams can also be defined for other types of sites, such as circles [66] and line segments [38]. A survey of Vornoi diagrams can be found in [7]. Figure 2.1: Voronoi diagram of a set of points. If the location of each point is not known precisely, but lies within a region of uncertainty, then we cannot construct the Voronoi diagram of the points. Instead, we define the guaranteed Voronoi diagram for the uncertain regions to be a set of cells, where the cell for site i consists of those points that are guaranteed, for each feasible set, to lie within the ith cell of the standard Voronoi diagram of the set. We denote the guaranteed Voronoi diagram of uncertain regions D by Ve (D). Figure 2.2 shows the guaranteed Voronoi diagram for a set of uncertain discs. One characteristic of guaranteed Voronoi diagrams that is immediately clear from Figure 2.2 is that for some points, we cannot guarantee a closest site. These points form a subset of the plane that we call the ‘neutral zone’. The Delaunay triangulation (or Delaunay graph) of a set of points is the dual of the Voronoi diagram of the points: an edge exists between two sites if and only if the Voronoi cells of the sites share an edge. The Delaunay triangulation has the property that for any triangle, the circumcircle of the triangle’s points contains no other points from the set in its interior (Figure 2.3). It has been shown that the Delaunay triangulation is the triangulation that maximizes the minimum angle of 1 Throughout this thesis, when discussing optimality, we will assume the standard algebraic decision tree model of computation. 12 Figure 2.2: Guaranteed Voronoi diagram for uncertain discs. its triangles. For this reason, they are usually the preferred triangulation of a set of points in graphics applications. Figure 2.3: Delaunay triangulation of a set of points, the dual of the Voronoi diagram (circumcircle for one triangle is shown). A set of imprecise points is associated with many feasible sets, and an edge that appears in the Delaunay triangulation of one feasible set may be absent in that of another feasible set. We define a guaranteed Delaunay edge to exist for a pair of uncertain regions drawn from a set D iff an edge associated with this pair exists 13 in the Delaunay triangulation of every feasible set of D. A Gabriel graph [27] is a subgraph of the Delaunay triangulation. A Gabriel edge exists between two points if the smallest disc that contains the points contains no other points in the set. We define a guaranteed Gabriel edge to exist for a pair of uncertain regions from D iff the corresponding edge exists in the Gabriel graph of every feasible set of D. 2.1.1 Related Work The guaranteed Voronoi diagram of a set of uncertain regions is closely related to the standard Voronoi diagram of those regions. Thus our results rely heavily on properties of standard Voronoi diagrams such as diagrams for discs ([66],[62]) and diagrams for segments [38]. There are two other generalizations of Voronoi diagrams that are closely related to our research. In an order-k Voronoi diagram, each cell is associated with a subset of the sites, and points within the cell have the property of being closer to each of the subset’s sites than to any of the sites not in the subset. These diagrams were first introduced by Miles [51], and were investigated from an algorithmic perspective by Shamos and Hoey [65]. In an additively weighted Voronoi diagram, each cell is associated with a single site, but the distance from a site to a point includes an additive weighting factor. The first algorithm to construct such diagrams was given by Drysdale and Lee [20], and improved algorithms were later given by Fortune [26] and Rosenberger [62]. If the weights are negative, then these diagrams are equivalent to the standard Voronoi diagram of discs, where each site is the origin of a disc with radius equal to the negative of the site’s weight.2 Voronoi diagrams involving uncertain sites that have a uniform probability distribution within discs were investigated in [8]. The authors defined the expected closest disc to a query point q to be the disc whose site has the smallest expected distance to q, and the probably closest disc to q to be the disc whose site is probably closest to q. These two definitions are not equivalent, and the authors have shown that while the Voronoi diagram for expected closest discs can be constructed us2 The equivalence holds as long as we allow distances to be negative, e.g. when a point is in a disc’s interior. 14 ing standard Voronoi algorithms, that of probably closest discs is more complex and cannot. Voronoi diagrams of probably closest discs have also been investigated in [4], and are shown to be useful for answering nearest neighbor queries in uncertain databases. Khanban investigates guaranteed Voronoi diagrams of uncertain convex polygons (with a small number of edges) in [36]. He points out that the cells of these diagrams are disjoint, and that the diagram contains a neutral zone. He proves that the diagrams are computable, but does not provide any complexity bounds or algorithms. Zone diagrams [5], like guaranteed Voronoi diagrams, have a neutral zone. In zone diagrams, for a point to be in a site’s region, it must be closer to the site than to any point in any other site’s region. The recursive nature of this definition raises the question of the uniqueness and existence of zone diagrams; a question that Asano et al [5] answered (positively). The fundamental operation in many algorithms for calculating Delaunay triangulations is an in-circle test, which determines whether a (query) point is inside or outside the circumcircle associated with a triangle’s points. Ely and Leclerc [24] investigated the application of interval arithmetic to this in-circle test for the case when the points are imprecise. They found that as the uncertainty in the input increases, the likelihood that the in-circle test will give reliable results decreases dramatically. To remedy this, they associated with each point a disc-shaped region of uncertainty, and showed that the set of all possible circumcircles consistent with the triangle points can be described by eight circles tangent to the triangle points’ discs, and that if the query disc lies outside (resp., inside) these circles, then the in-circle test can reliably report false (resp., true). They found that an in-circle test that uses this tangent method is about sixteen times less likely to fail to report an answer than one that only uses interval arithmetic. Khanban and Edalat [37] have argued that in practice, imprecision in input points are best represented by rectangles rather than discs, since the points are most often reported as a pair of x and y coordinates. They show how the uncertain rectangles corresponding to the three vertices of each Delaunay triangle give rise to six tangent circles, which can be used for an in-circle test in a similar manner to that of Ely and Leclerc. Khanban gives an algorithm for constructing the guaranteed 15 Delaunay edges of uncertain convex polygons (with a small number of edges) that that runs in expected O(n log n) time [36]. In a different take on imprecision, Löffler and Snoeyink [46] investigated the problem of preprocessing a set of disjoint uncertain unit discs to later allow efficient construction of the Delaunay triangulation of a feasible set of the discs. They show that once O(n log n) time has been spent preprocessing the uncertain discs, the Delaunay triangulation of any feasible set from the discs can be generated in O(n) additional time. Randomized algorithms that extend and simplify these results are given by Buchin et al in [13]. For the case of uncertain discs with the property that no point lies in more than k discs, they provide an algorithm that (after an O(n log n) preprocessing step) can generate the Delaunay triangulation of any feasible set in expected O(n log k) time. The algorithms in both of these works are largely of theoretical interest, since the constants hidden in their postprocessing step running times are large enough that it is often faster just to spend O(n log n) time constructing the Delaunay triangulation of each feasible set using a standard algorithm. By contrast, Devillers [19] has shown that for disjoint uncertain discs of unit radius, there exists a much simpler randomized algorithm that produces better (expected) postprocessing running times than the standard Delaunay triangulation algorithms. The tolerance of a geometric or combinatoric object is a concept that was introduced by Abellanas et al [1]. It refers to the object’s ability to withstand perturbation while still maintaining some property. For example, the tolerance of the Delaunay triangulation of a set of points S is the maximum amount by which the points can move and still have the same Delaunay triangulation. It was shown in [2] that determining for this problem can be done in O(n) time (if the Delaunay triangulation is given), or in Θ(n log n) time (otherwise). This problem is related to our work in the following way: the tolerance of the Delaunay triangulation of S is the maximum radius that each of a set of uncertain discs D centered at S can have and still satisfy the property that every Delaunay edge of S is a guaranteed Delaunay edge of D. 16 2.1.2 Contributions This chapter is concerned with Voronoi diagrams, Delaunay triangulations, and Gabriel graphs of imprecise points.3 We prove tight bounds on the complexities of the guaranteed Voronoi diagrams of uncertain discs and uncertain polygons, and present worst-case optimal algorithms for their construction. We define a new diagram, the possible Voronoi diagram. It reports the smallest subset of regions that are guaranteed to contain the closest site to a query point. In addition to proving a tight bound on its complexity, we provide a construction algorithm that is worst-case optimal up to logarithmic terms. Finally, we provide worst-case optimal algorithms for determining both the guaranteed Delaunay edges and the guaranteed Gabriel edges associated with a set of uncertain discs. 2.2 Guaranteed Voronoi Diagram In this section, we will investigate some properties of guaranteed Voronoi diagrams. We will then define a useful mapping function between boundary points of standard and guaranteed Voronoi cells. We will then describe algorithms for constructing these diagrams for uncertain discs and (disjoint) uncertain polygons. 2.2.1 Properties We denote the standard Voronoi diagram of a set of regions D by V (D), and the cell corresponding to region i by Ri . We denote the order-k Voronoi diagram of a set of regions D by V k (D), and a cell of this diagram corresponding to the |S| = k regions {Di∈S⊆{1...n} } by RS . Hence V (D) = V 1 (D), and Ri = R{i} . Suppose we have a set of uncertain regions D. The bisector of a pair of regions Di and Dj is the set of points that are equidistant to both regions. We denote this bisector by hhi, jii. Formally, hhi, jii = {p | min d(p, x) = min d(p, y)} x∈Di y∈Dj = {p | d(p, Di ) = d(p, Dj )} . 3 Portions of this chapter appeared in [63]. 17 Let H(i, j) be the set of points that are guaranteed to be at least as close to site i as site j. That is, H(i, j) = {p | ∀x ∈ Di ∀y ∈ Dj d(p, x) ≤ d(p, y)} . We denote the boundary of H(i, j) by hi, ji; formally, hi, ji = {p | max d(p, x) = min d(p, y)} . x∈Di y∈Dj ei , is: The guaranteed cell for site i, denoted R ei = R \ H(i, j) . (2.1) j6=i ei form the guaranteed Voronoi diagram for the The boundaries of all such cells R set D, and we denote it by Ve (D). An edge of hi, ji in Ve (D) is a maximal connected ei . set of points p ∈ hi, ji that lie on the boundary of cell R In order to further investigate the properties of Ve (D), we will require some additional terminology. Let A be a disc. Shrinking A refers to the process of reducing the radius of A while keeping its origin fixed. If p is a point on A’s boundary, then shrinking A with respect to p refers to the process of reducing the radius of A while simultaneously moving its origin towards p, so that p remains on the boundary of (the shrinking) A. If p is a tangency point of A and some region B, then shrinking A with respect to p ensures that (the shrinking) A and B will remain tangent at p. When clear from the context, we refer to such an operation as shrinking A with respect to B. ei in Ve (D) iff there exists a disc Cp centered at p such Lemma 2.1 Point p is in R that (i) Cp has inside-tangent Di , and (ii) D \ {Di } does not penetrate Cp . In ei iff Cp has outside-tangent some Dj6=i ∈ D. addition, p is on the boundary of R ei and H(i, j). Proof. This follows immediately from the definitions of R Figure 2.4 illustrates Lemma 2.1. e i = Ri . Lemma 2.2 If region Di of D is a single point, then R 18 p Di ei . Figure 2.4: Lemma 2.1: p ∈ R Proof. Suppose Di is a single point. We can apply an analog of Lemma 2.1 (one concerned with standard Voronoi diagram cells) to show that for every point p ∈ Ri , there exists a disc Cp that has outside-tangent Di , and is not penetrated by D \ {Di }. Since Di is a point, each such Cp also has inside-tangent Di , satisfying Lemma 2.1. ei of Ve (D) is a subset of the corresponding cell Ri in the Lemma 2.3 Every cell R standard Voronoi diagram V (D). ei , then the distance from p to the farthest point of Proof. If p is a point within R Di is not greater than the distance from p to the nearest point of any region D6=i . Since the former distance is an upper bound on the distance from p to the nearest point of Di , p must lie within Ri as well. ei and R ej will Lemma 2.4 If Di and Dj are intersecting regions of D, then both R be empty. Proof. Let Di and Dj be intersecting regions of D, and p be any point in the plane (Figure 2.5). Observe that there must exist a feasible set of D in which si is closer ej is empty. A symmetric argument implies that to p than is sj , which implies that R ei is empty as well. R A region whose cell is empty can still influence the cell of another region. For example, the two uncertain discs in the bottom right of Figure 2.2 have empty cells, yet still induce edges in some nearby cells. 19 p si sj Di sj si Dj Figure 2.5: Lemma 2.4. ei of Ve (D) is convex. Lemma 2.5 Each cell R ei . We show that segment pq lies within Proof. Let p and q be arbitrary points in R ei . Lemma 2.1 implies that there exist discs Cp centered at p and Cq centered at q R that contain Di , and are penetrated by no disc D6=i . This implies that Di ⊆ Cp ∩Cq . Now, each point p0 ∈ pq is the center of a disc Cp0 that contains Cp ∩ Cq , and lies within Cp ∪ Cq ; hence Cp0 satisfies condition (ii) of Lemma 2.1. We can shrink Cp0 until it has inside-tangent Di , at which point it will satisfy condition (i) as well. ei . Thus p0 ∈ R 2.2.2 A Mapping between Boundaries of Voronoi Cells We will make extensive use of a mapping from boundaries of cells of Ve (D) (i.e., the guaranteed diagram) to boundaries of cells of V (D) (i.e., the standard diagram). We will construct the mapping, δ(·), in the following way. Consider a point p on ei in Ve (D), and its disc Cp from Lemma 2.1. Since p is on the the boundary of R boundary of the region, Cp has outside-tangent some Dm6=i at some point a;4 see Figure 2.6. We now shrink Cp with respect to a, until it reaches Cp0 that has outside-tangent 4 If more than one region Dm is outside-tangent to Cp (or if Dm is tangent to Cp at more than one point), then we make δ(p) well-defined by selecting point a according to some total order on possible a’s. 20 Dm a p0 p Cp0 Di Cp Figure 2.6: Definition of δ(p). Di . Now, p0 is on an edge of hhi, mii on the boundary of Ri ⊆ V (D), and we define δ(p) to be this point p0 . If we start at a point s on the boundary of a connected region, and traverse the boundary in a ccw direction, then we denote a point p being encountered before q during this traversal by p ≺s q. ei within Ve (D), and Lemma 2.6 If p, q, and s are points on the boundary of cell R p ≺s q, then δ(p) δ(s) δ(q). Proof. We will first show that if the δ(·) mapping does not preserve this ordering, then some pair of the line segments pδ(p), qδ(q), sδ(s) must intersect; we will then derive a contradiction from this fact. ei , and L be segment pδ(p) (Figure 2.7). Now, L Let p be a boundary point of R ei ; for if it did, some point p0 interior to L must lie in the interior does not penetrate R ei , but no disc centered at p0 can contain Di without also being penetrated by of R some D6=i . Note also that L is within Ri , since for each p0 ∈ L, (i) there exists a disc centered at p0 that intersects Di and that has outside-tangent some D6=i , and (ii) no smaller disc at p0 intersects any D6=i . 21 Di Ri L f R i p δ(p) Figure 2.7: Lemma 2.6. Thus if p ≺s q and δ(p) δ(s) δ(q), then p 6= q, δ(p) 6= δ(q), and some two of the segments sδ(s), pδ(p), qδ(q) intersect. Without loss of generality, assume pδ(p) intersects qδ(q) at a point v. By Lemma 2.1, there exists disc Cp which has outside-tangent some D6=i at a, such that δ(p) ∈ pa. Similarly, disc Cq exists which has outside-tangent some D6=i at b(6= a) where δ(q) ∈ qb. The mapping δ(·) implies that every point on pδ(p) is the center of a disc that is penetrated by no D6=i , and has outside-tangent some D6=i at a. Since a symmetric situation holds for q, δ(q), and b, v must be the center of a disc Cv that has a and b on its boundary. Now observe that every point of Cv except a must lie in the interior of Cp , including b (Figure 2.8); hence Cp is penetrated by some D6=i , a contradiction. 2.2.3 Uncertain Discs We now consider the case where the uncertain regions are discs (e.g., Figure 2.2). Each disc has a nonnegative radius ri , and a center oi . Every point p ∈ hi, ji 22 D6=i Cp a Cv pδ(p) v p qδ(q) b D6=i Figure 2.8: Lemma 2.6: pδ(p) and qδ(q) cannot intersect. satisfies d(p, oi ) + ri = d(p, oj ) − rj . (2.2) Since ri and rj are constants, the points p which satisfy (2.2) lie on the arm of the hyperbola with foci at oi and oj that is closest to oi (if Di and Dj intersect, then hi, ji = ∅). ei within Ve (D) may contain two or more distinct edges The boundary of a cell R of hi, ji; see Figure 2.9. We now show that the number of edges in the guaranteed Voronoi diagram of a set of uncertain discs is at most twice that of the standard Voronoi diagram of the discs. We do this by showing that each edge of each cell of the guaranteed diagram maps to a distinct edge in the corresponding cell of the standard diagram. The doubling factor is explained by the fact that each edge of the standard diagram is shared by two cells, whereas edges of the guaranteed diagram typically are not. Theorem 2.7 The guaranteed Voronoi diagram of n uncertain discs has O(n) edges. Proof. We will show that the number of edges in Ve (D) is at most twice the number of edges in V (D). The theorem then follows from the fact that V (D) has O(n) 23 Dj Dk hi, ki hi, ji hi, ji Di f R i ei . Figure 2.9: Multiple edges of hi, ji on the boundary of R edges [66]. ei in ccw order. Lemma 2.5 ensures that these Consider the edges around R edges are adjacent. We charge each edge E of hi, ji ∈ Ve (D) to the edge F of hhi, jii on which δ(p) lies, for p the ccw-first point of E (or any interior point p if E is ccw-infinite). Suppose two distinct edges E1 and E2 of hi, ji map to the same edge F of hhi, jii in Ri ⊆ V (D). Since E1 and E2 are distinct but both of hi, ji, there must exist an edge E 0 of hi, ki (with k 6= j) between them in the ccw ei that maps to some other edge F 0 of hhi, kii in Ri . This contradicts traversal of R Lemma 2.6 since all points of F either precede or follow the points of F 0 in ccw order. Thus each edge in V (D) of hhi, jii is charged at most twice: once by an edge of hi, ji, and once by an edge of hj, ii. Hence Ve (D), like V (D), has O(n) edges. We now show how Ve (D) for a set of discs can be constructed by first constructing V (D) for the discs, then performing a linear-time transformation from V (D) to Ve (D). Let Iei be the sequence of pairs (i, j) corresponding to the edges hi, ji encoun24 ei , when starting from an edge contered in a ccw traversal of the boundary of R taining some point s on the boundary. We similary define Ii to be the sequence of pairs (i, j) corresponding to the edges hhi, jii traversed on the boundary of Ri , when starting from the edge containing δ(s). ei ∈ Ve (D), Iei is a subsequence of Ii . Lemma 2.8 For every cell R Proof. If (i, j) is in Iei , then since δ(·) maps points on edges of hi, ji to points on edges of hhi, jii, (i, j) must be in Ii . Furthermore, the order of pairs in Iei is preserved in Ii since δ(·) preserves this order by Lemma 2.6. One complication arises if the discs are not partially disjoint. In order to discuss this situation, we require some additional terminology. A disc Di ∈ D is an outer disc if it is not contained by any other disc of D. In addition, if Di is the largest (outer) disc containing Dj , then we will say that Dj is the child of Di , and that Di is its parent. Note that while a child may be contained by several discs, it has a single parent, and is a child of only that parent.5 Each child disc will have an empty region within V (D), since every point in the plane will be closer (allowing for negative distances) to the child’s parent than it is to the child. The standard Voronoi diagram thus offers no evidence that one disc contains another, which may lead to problems since (by Lemma 2.4) the parent’s guaranteed cell is actually empty. With these facts in mind, we are ready to describe our algorithm for constructing the guaranteed Voronoi diagram of uncertain discs. Theorem 2.9 Ve (D) for n uncertain discs can be constructed in O(n log n) time. Proof. We construct V (D) for the disc sites in O(n log n) time, using Fortune’s algorithm [26]. That algorithm can determine, as a side effect, (i) the set of outer discs, (ii) which outer discs are parents, and (iii) the set of children of each parent. By Lemma 2.4, every parent will have an empty cell within Ve (D), as will each of its children. Thus, we need construct guaranteed cells only for outer discs that are not parents. 5 In the case where two discs of equal size contain a particular disc, we will arbitrarily choose one of the containing discs to be its parent. 25 For each such disc Di , we construct its guaranteed cell from its standard cell Ri as follows. We generate the sequence Ii of sites comprising the boundary of ei Ri in linear time by a simple traversal. We then construct the boundary of R by generating and intersecting the sequence of hyperbolic arcs hi, ji for each pair (i, j) ∈ Ii . If any of these hyperbolic arcs do not exist, then Di and Dj must ei is empty. Lemma 2.8 ensures that intersect, which implies (by Lemma 2.4) that R ei ; hence, by we consider a correctly ordered super-sequence of the arcs bounding R ei can be constructed from performing a straightforward clipping procedure, cell R this super-sequence in time proportional to the length of Ii . Observe that if Dj is the parent of Dk , and Di contains Dk as well, then our algorithm will still attempt to construct a guaranteed cell for Di . It will correctly produce an empty guaranteed cell, since Dj must intersect Di (and hence hi, ji will not exist), and Dj will be a neighbor of Di within the standard Voronoi diagram.6 Since each of the O(n) edges of V (D) appears in two cell boundaries, the running time for the construction of the edges of all cells of Ve (D) is O(n). The time to construct Ve (D) is thus dominated by the time to construct V (D). The algorithm of Theorem 2.9 is worst-case optimal. To see this, observe that when the disc radii are all zero, Ve (D) is the standard Voronoi diagram of n points, which requires Ω(n log n) time to construct [65]. 2.2.4 Uncertain Polygons In this section, we present an efficient algorithm for constructing the guaranteed Voronoi diagram of a set of disjoint7 uncertain polygons (Figure 2.10). As noted earlier, this is an important class of uncertain region, since it includes axis-aligned squares, a type of region likely to be encountered in practice. The algorithm will take the same approach as that for the case of uncertain discs: namely, it will start by constructing V (D), the standard Voronoi diagram of the polygons, then ei . Our will transform each cell Ri of this diagram to its guaranteed counterpart, R analysis of the algorithm’s running time will also prove that these diagrams have linear complexity. 6 Actually, Dj might not be a neighbor of Di , but only if some other (outer) disc that intersects 26 Figure 2.10: Guaranteed Voronoi diagram of uncertain polygons. Throughout this section, we assume that D is a set of disjoint uncertain polygons that satisfies the following general position conditions: no vertex is contained in another polygon’s edge, no three edges intersect at a common point, and no disc exists that contains four points that are either polygon vertices, or points where the disc is tangent to a polygon edge. To further simplify the following discussion, we will also assume that the cells ei are both bounded. The algorithm can be easily modified to deal with Ri and R unbounded regions, without affecting its running time. Alternatively, we could view D as being a cluster of polygons centered within a large rectangle, and add a small triangle to D at each corner of this rectangle. This would ensure that the Voronoi cells of the original polygons are bounded, but would not change the structure of these cells except near the box boundaries. The boundary of each cell Ri is composed of subsets of bisectors hhei , ej ii, where ei is an edge of polygon Di , and ej is an edge of a polygon Dj6=i ∈ D. Di is. 7 Section 2.6 explains why we require the polygons to be disjoint. 27 To aid in our subsequent analysis, we note that each boundary point p on such a bisector is the center of a disc Cp that is not penetrated by any polygon of D, and has outside-tangent8 edges ei and ej of polygons Di and Dj respectively. The boundary of a guaranteed Voronoi cell is similarly composed of subsets of (guaranteed) bisectors hv, ei, where v is a vertex of Di , and e is an edge of some Dj6=i ∈ D. Each point p on such a bisector is the center of a disc Cp which has inside-tangent Di at v, has outside-tangent Dj at e, and is not penetrated by any ei . polygon D6=i ∈ D. We will refer to Cp as a guaranteed disc for R ei can As a consequence of our general position assumption, each vertex p of R be identified as a triple hv1 , v2 , ei or a triple hv, e1 , e2 i. The former corresponds to a guaranteed disc Cp that has inside-tangent Di at vertices v1 and v2 , and has outside-tangent edge e, while the latter corresponds to a guaranteed disc Cp that has inside-tangent Di at v, and has outside-tangent distinct edges e1 and e2 . In both cases, we order the triples to agree with the ccw order of the tangent points of these vertices and edges with Cp . ei that lies on bisector hv, ei will terminate at a vertex Lemma 2.10 Each edge of R of the form hv, v 0 , ei or hv, e, e0 i. ei is the center of a Proof. Each point p ∈ hv, ei that lies on the boundary of R guaranteed disc Cp in a family of such discs tangent to v and e. Observe that the arc on the boundary of each Cp from v ccw to (the point of tangency with) e is receding, in the sense that the interior points of this arc will not appear in subsequent discs in the family. Similarly, the boundary arc from e ccw to v is advancing, in the sense that points interior to the arc will be in the interior of the subsequent discs in the family (unless they appear on the receding arc of a later disc, which may occur as the point of tangency with e changes). Now, note that for ei , some vertex of Di must appear on the receding arc of Cp , or p to be a vertex of R some edge of a polygon D6=i ∈ D must appear on the advancing arc of Cp . Let Vi represent the ccw-ordered sequence of vertices on the convex hull of polygon Di , and Ei represent the ccw-ordered sequence of edges e that appear in 8 We can consider tangency with an edge endpoint as a special case of tangency with the edge itself, and to eliminate ambiguity, we can define each (directed) edge as closed at its source endpoint, and open at its destination endpoint. 28 the bisectors hh·, eii that form the boundary of cell Ri . We will let k denote the total complexity of Vi and Ei . Lemma 2.11 The sequence of vertices of uncertain polygon Di that induce verei are a subsequence of Vi . tices of R Proof. We will denote point u preceding point v on a (ccw oriented) region boundei , with p ≺ q ≺ r. ary by u ≺ v. Let p, q, and r be distinct boundary points of R ei is convex, r must lie to the left of pq. Since R ei , there exists a guaranteed disc Cs that For each point s on the boundary of R has inside-tangent Di at a vertex ts . Clearly, ts must be on the hull of Di ; hence, we need only show that the ordering of these vertices agrees with that of the vertices ei . of R Now, Di must lie in the common intersection I = Cp ∩ Cq ∩ Cr . We can show (by our general position assumption, and the fact that each disc intersects, but is not penetrated by, some D6=i ) that none of the three discs is contained by the union of the other two. This implies that the boundary of I is comprised of a single circular arc from each of Cp , Cq , and Cr , and that these arcs have a ccw ordering consistent with that of p, q, and r (Figure 2.11). Since each disc’s tangent point must lie on its corresponding arc, this implies that tp tq tr ; and as this holds for any choice ei , the proof is complete. of distinct points p ≺ q ≺ r on the boundary of R Lemma 2.12 If disc Cp contains polygon Di and is penetrated by any edge of D \ {Di }, then it will be penetrated by some edge of Ei . Proof. Suppose Di ⊂ Cp , and Cp is penetrated by edge e of U = D \ {Di }. Starting with C = Cp , we can shrink it until it has outside-tangent either Di or U , then shrink it further with respect to that tangency point until it has outside-tangent both Di and some edge e of U , and is not penetrated by D. Now, the final C is a guaranteed disc, hence e ∈ Ei ; and since it lies within Cp , e penetrates Cp . ei in two stages. The first stage finds a guaranteed disc for We will construct R ei . If none exists, then R ei must be empty. Otherwise, its origin is a point on the R ei , and starting from this point, the second stage generates the ordered boundary of R ei . sequence of vertices of R 29 tp Ri tr f r Di p Cr tq q I Cp Cq Figure 2.11: Lemma 2.11. We begin the first stage by constructing the smallest disc containing polygon Di . If this disc is not penetrated by any other polygons of D, then we perform the following grow operation: we expand the disc, maintaining its inside-tangency with Di , until it has outside-tangent some D6=i , at which time it will be a guaranteed disc. If instead, the disc is penetrated by some edge of a polygon D6=i , then we perform a move operation: we transform the disc, maintaining its tangency with Di , until no edges penetrate the disc (or we determine that no guaranteed disc exists). We now describe these steps with more detail. We first construct the farthest vertex Voronoi diagram of Vi , and the nearest edge Voronoi diagram of Ei , in O(k log k) time. Using the farthest vertex Voronoi diagram, we can find Cs , the smallest bounding disc of Vi , in O(k) time [65]. We query the nearest edge Voronoi diagram to determine the nearest edge e ∈ Ei to s. If e does not penetrate Cs , ei . then by Lemma 2.12, no edge of D \ {Di } will either, which implies that s ∈ R To perform the grow operation, we start with any vertex v ∈ Vi that lies on the boundary of Cs . We then determine (in linear time, by iterating over Ei ) the small- 30 → est disc Cp that is (i) centered at a point on ray − vs, (ii) tangent to v, and (iii) tangent ei , by Lemma 2.12. to some edge of Ei . Now, Cp must be a guaranteed disc for R v Cs Di s u e w ei . Also shown are farthest point Figure 2.12: Finding point on boundary of R Voronoi diagram of Vi (solid lines) and bisector of u and w (dashed line). If the bounding disc Cs is penetrated by edge e of some other polygon, we perform the move operation instead. Its steps are a little more involved, and will require the following lemma. Lemma 2.13 ([18], Lemma 9.4) If disc C with boundary points α and β intersects − → points p and q lying on opposite sides of directed line αβ, then every other disc intersecting α and β will be penetrated by either p or q (Figure 2.13). Let u and v be two of the three vertices of Vi that lie on the boundary of Cs , such → (Figure 2.12). that some portion of e lies in the interior of C to the left of − uv s → maintaining its We transform Cs by moving its center towards the right of − uv, tangency with u and v. In a manner analogous to that of the proof of Lemma 2.10, the boundary of Cs can be partitioned into advancing and receding arcs, separated → We transform C until one of three events occurs (Figure 2.14): by the ray − uv. s 31 β C p q α Figure 2.13: Lemma 2.13. Event 1: e does not penetrate Cs . Just prior to this point, Cs was penetrated by some edge of Ei ; hence, it must now have outside-tangent one of these edges on its receding arc. Since it also has inside-tangent Di , we can return Cs as a guaranteed disc. Event 2: Some other vertex w ∈ Vi appears on the receding arc of Cs . We must replace either u or v with w before continuing to transform Cs , otherwise w will not be contained by Cs . If e penetrates the portion of Cs lying to the left of uw (as in Figure 2.14), then we replace v with w; otherwise, we replace u with w. Afterwards, we continue transforming Cs . → Event 3: e penetrates Cs to the right of − uv. → by the Observe that just prior to this point, Cs was penetrated to the left of − uv → by the new edge e (which may be the same previous edge e and to the right of − uv as the previous e, if the last event was of type 2). It follows (by Lemma 2.13 and ei is empty, so we return null. Lemma 2.1) that R Clearly, we will want to grow Ds in a finite number of discrete steps. At each iteration, it will suffice to construct two candidates for the next Cs , both tangent to u and v: a type 1 candidate, which has outside-tangent edge e; and a type 2 candidate, which is tangent an additional vertex w of Di . We replace Cs with the candidate whose origin is closest to the current s. Next, we update e to be the nearest edge of Ei to the new location of s. We then determine which of the three events applies to the new Cs , and take the appropriate action. 32 e Cs Di Cs Di s Di s s u e v u e v Cs u v w Event 1 Event 2 Event 3 Figure 2.14: Events that can occur while growing Cs . Let us examine the running time of an iteration of the above procedure. At the start of each iteration, we can determine e, the closest edge of Ei to s, by querying the nearest edge Voronoi diagram. This takes O(log k) time. We can determine which of the three events has occurred, and take the appropriate action, in constant time (note that an event of type 2 has occurred iff a type 2 candidate was chosen to replace Cs in the previous iteration). Constructing type 1 candidates can be done in constant time. To construct a type 2 candidate, observe that s is moving along an edge f of the farthest point Voronoi diagram bordering the cells of u and v, and that the next (Voronoi) vertex encountered will be the center of the type 2 candidate. If we maintain a pointer to f , then the candidate and its associated vertex w can be determined in constant time. The running time of each iteration is thus dominated by the time spent updating e, which can be done in O(log k) time. To bound the number of iterations, let us consider the sequence of discs Cs generated by this procedure. Observe that (i) each iteration sees a vertex of Vi or an edge of Ei becoming outside-tangent to the receding arc of some Cs , and (ii) the receding arc for a disc is a strict subset of that of the previous disc in the sequence. We can thus charge each iteration to a unique element of Vi or Ei . It follows that there are at most k iterations, and we can make the following claim. ei can Lemma 2.14 Given Ei and Vi of total complexity k, a boundary point of R ei is empty) in O(k log k) time. be found (or it can be determined that R 33 ei If we apply the above procedure and determine that no guaranteed disc for R ei must be empty. Otherwise, we are ready to perform the second exists, then R stage, which actually constructs the guaranteed cell for Di . We will make use of a circular array data structure, one that provides the constant time operations peek() to return the current item, and advance(), to advance the cursor to the next item (or to the first item, if the current item was the last item). We perform these steps: 1. Determine the bisector hv, ei that contains p, the origin of the guaranteed disc Cp found earlier. Also, determine which edge of Ei contains δ(p) (the corresponding boundary point of Ri ; see Section 2.2.2), by finding the disc Cδ(p) which is outside-tangent (at the same point) the same edge of Ei as Cp . 2. Construct (i) circular array AV , consisting of the ordered vertices Vi , with its cursor initially pointing to the vertex following v; and (ii) circular array AE , consisting of the ordered edges Ei ,9 with its cursor initially pointing to the edge following that which contains δ(p). 3. Construct candidate vertices hv, v 0 , ei and hv, e, e0 i, as follows. Let v 0 = AV .peek() and e0 = AE .peek(). If v = v 0 (resp., e = e0 ), advance AV (resp., AE ), and repeat from step 3. Otherwise, construct the discs hv, v 0 , ei and hv, e, e0 i. Let Cp0 be the disc of these two whose origin represents the shorter distance forward from p along hv, ei. 4. If Cp0 does not contain Di , then advance AV , and repeat from step 3. 5. If Cp0 is penetrated by any polygon D6=i ∈ D, then advance AE and repeat from step 3. 6. If p0 was already output as the first vertex, then stop. Otherwise, output p0 , replace p with p0 , replace v or e with the new vertex or edge that produced Cp0 , and repeat from step 3. Let us now examine the above steps in detail, to prove that the procedure genei . Point p lies on bisector hv, ei, and serves as our starting point for genererates R ei . ating the vertices of R 9 Consecutive occurrences of the same edge can be omitted. 34 ei is either of the form hv, v 0 , ei or hv, e, e0 i. By Lemma 2.10, the next vertex of R In the former case, Lemma 2.11 implies that v 0 is either the next element of Vi , or appears after some number of elements of Vi that can be discarded. In the latter case, we can use an argument analogous to that of the proof of Lemma 2.8 to show that if e0 is not the next element of Ei , then it must appear after one or more elements of Ei (which can similarly be discarded). Once we have chosen a candidate next vertex p0 in step 3, we verify that its associated guaranteed disc Cp0 contains Di . If not, then some vertex w on the hull of Di lies outside of Cp0 . Observe that (i) w must have appeared on the receding arc of some disc Cq that is tangent to v, w, and e; and (ii) v 0 must lie within Cq and (since v 0 is the next vertex on the convex hull of Di ) to the right of vw (Figure 2.15). ei must be the center of a disc that is not penetrated Note also that every point in R by e, but contains both v and w. This implies that v 0 cannot lie on the boundary of such a disc. Hence, we can discard v 0 (by advancing AV ), and repeat from step 3. C p0 Di e v Dj v0 Cq w ei (Cp0 equals hv, v 0 , ei in this Figure 2.15: v 0 cannot contribute to R example). We next verify that no other polygon penetrates Cp0 in step 5. If some peneei . It follows that there exists a point q trating edge f exists, then p0 cannot lie in R 35 between p and p0 on hv, ei such that f lies on the advancing arc of Cq . This implies ei . Hence, by Lemma 2.6, we can discard e0 (by that hv, e, f i is the next vertex of R advancing AE ) and repeat from step 3 (actually, we can safely advance AE until f appears). By the time we reach step 6, we have established that Cp0 is a guaranteed disc ei . Observe that the next vertex of ei , which implies that p0 is a vertex of R for R ei cannot be hv, w, ei where w is a vertex of Vi other than v 0 , since each such R w must appear after v but before v 0 on the convex hull of Di (by Lemma 2.11), and will thus have been eliminated as a possibility by some execution of step 4. ei cannot be Similar reasoning can be applied to show that the next vertex of R hv, e, f i where f is an edge of Ei other than e0 , for each such edge f will have ei . If p0 been eliminated in step 5. Hence, p0 must in fact be the next vertex of R was already output as the first vertex, then we have generated the entire boundary ei , and we are done. of R Now that we have shown that the algorithm correctly generates Ve (D), we will determine its running time. In order to convert standard Voronoi cells to their guaranteed counterparts, we must construct the standard Voronoi diagram of the polygons. This can be done in O(n log n) time (e.g., [26, 38, 76]). The first ei , can be done in O(k log k) time, stage, finding an initial guaranteed disc for R by Lemma 2.14. Note that the algorithm we gave as proof of this lemma constructs two Voronoi diagrams: the farthest vertex of Vi , and the nearest edge of Ei . We will make use of these diagrams below. Steps 1 and 2 run in O(k) time, and step 3 runs in constant time. We can test if Cp0 contains Di by querying the farthest vertex Voronoi diagram to verify that v is the farthest point of Di from p0 . Hence, step 4 runs in O(log k) time. In step 5, we need to determine if any D6=i penetrates Cp0 . We can do this by querying the nearest edge Voronoi diagram to see if the nearest edge of Ei penetrates Cp0 . If not, then by Lemma 2.12, no polygon of D \ {Di } penetrates C. Step 5 thus runs in O(log k) time. Since step 6 runs in constant time, a complete iteration of steps 3 through 6 runs in O(log k) time. Note that each such iteration either (i) eliminates a vertex of ei . Hence, at AV , (ii) eliminates an edge of AE , or (iii) outputs a unique vertex of R most O(k) iterations are performed, for a total running time of O(k log k) per cell. 36 Since the total of k for all the guaranteed cells is O(n), we can claim that: Theorem 2.15 The guaranteed Voronoi diagram of disjoint uncertain polygons with a total of n vertices has O(n) complexity, and can be constructed in O(n log n) time. 2.3 Possible Voronoi Diagram One of the significant differences between guaranteed Voronoi diagrams and their standard counterparts is the presence of a neutral zone: areas of the plane for which no guaranteed closest site can be determined. In this section, we introduce a Voronoi diagram of imprecise points which has no neutral zone, and thus can make a guarantee about every point in the plane. We achieve this by modifying the guarantee. In particular, this variant determines, for each point p, the set of possibly closest sites: the smallest set of sites such that one is guaranteed to be a closest site to p (throughout this section, we will assume that the uncertain regions are discs). We call the resulting partition of the plane, which we denote by Ve {} (D), a possible Voronoi diagram (Figure 2.16). In contrast to guaranteed Voronoi diagrams, these diagrams have cells that are not necessarily connected. In Figure 2.16, for instance, the two shaded regions belong to the same cell. The possible Voronoi diagram can also be defined by generalizing equation (2.1): eS = R [h \ i [ eS 0 H(i, j) − R S 0 ⊂S i∈S j ∈S / e∅ = ∅. Note that only the closure of this definition of R eS is strictly where R equivalent to the first definition. The following construct will prove useful in generating Ve {} (D). The possible cell for a site i, denoted Pei , is the set of all points for which it is possible that site i is the closest site. Formally, Pei = {p | ∃q1 ∈ D1 , q2 ∈ D2 , . . . , qn ∈ Dn ∀j d(p, qi ) ≤ d(p, qj )} . The possible cell for one particular site is shown in Figure 2.17. 37 Figure 2.16: Possible Voronoi diagram. Pfi Figure 2.17: The possible cell for a site. The possible cells for a set of sites do not generally have disjoint interiors, unless each uncertain region is a single point, in which case the possible cells are just standard Voronoi cells for point sites (and each site’s location is no longer uncertain). Lemma 2.16 Point p is within Pei iff there exists a disc Cp centered at p that intersects Di , and no D6=i ∈ D is strictly interior to the disc. In addition, p is on the boundary of Pei iff Cp has outside-tangent Di and inside-tangent some D6=i . 38 Observe that the possible Voronoi diagram associates each point in the plane with the subset of sites whose possible cells contain the point. Hence, to construct the possible Voronoi diagram, it suffices to construct the arrangement of the possible cell boundaries. To do this efficiently, we will make use of the additively weighted Voronoi diagram of a set of point sites derived from the set of discs. We associate with each disc Di a weighted site (oi , −ri ). We denote the resulting Voronoi diagram by V m (D), and the cell for (oi , −ri ) by Rim . The bisector between sites i and j is hhi, jii = {p | d(p, oi ) + ri = d(p, oj ) + rj } To gain some intuition concerning this type of Voronoi diagram, note that it is essentially a standard Voronoi diagram of discs, except that the distance from a point to a disc is defined to be the distance from the point to the disc’s most distant boundary point. Let us replace a particular site i’s weight with +ri . Any point p on edge hhi, jii bordering cell i now satisfies d(p, oi ) − ri = d(p, oj ) + rj . By Lemma 2.16, this implies that the set of all edges bordering the new cell i forms the boundary of Pei . In addition, since the new diagram is still an additively weighted Voronoi diagram, each possible cell has linear complexity. Observe that the possible cell is produced by modifying a single site’s additive weight in such a way that the possible cell is a superset of the corresponding cell of V m (D). It follows that the possible cell’s vertices must lie on existing edges of V m (D); see Figure 2.18. This suggests an approach to generating the possible cells: we construct the additively weighted Voronoi diagram of sites {(oi , −ri ) | i = 1 . . . n}, then for each disc Di in turn, we: 1. insert a new site (oi , +ri ), 2. extract Pei from this new site’s cell boundary, and 39 Di Pfi Figure 2.18: Possible cell vertices lie on edges of V m (D). 3. delete (oi , +ri ). As we will see, we will not actually perform this insertion and deletion of sites, and can instead extract each Pei by performing a search along edges of V m (D). Before describing the algorithm in more detail, we need to investigate further the relationship between intersecting discs, cells of V m (D), and possible regions. Lemma 2.17 If Di intersects Dj , then Rjm ⊆ Pei . Proof. If Di intersects Dj , then for any point p ∈ Rjm , we have, for all Dk , d(p, ok ) + rk ≥ d(p, oj ) + rj ≥ d(p, oi ) − d(oi , oj ) + rj (triangle inequality) ≥ d(p, oi ) − (ri + rj ) + rj = d(p, oi ) − ri , which implies that p ∈ Pei . Observe that Lemma 2.17 implies that for every Di ∈ D, Rim ⊆ Pei . 40 Lemma 2.18 If Di contains Dj , then Rim is empty. Proof. If Di contains Dj , then the most distant boundary point of Di from any point p will be farther from p than the most distant boundary point of Dj . Lemma 2.19 If Di and Dj are nonintersecting discs, and Rjm is nonempty, then Rm 6⊆ Pei . j Proof. Suppose Rjm is nonempty. Lemma 2.18 implies that oj is closer to the most distant point of Dj than it is to the most distant point of any D6=j ; hence oj ∈ Rjm . Since Di does not intersect Dj , oj is closer to the most distant point of Dj than it is to any point of Di ; hence oj ∈ / Pei . Our first step is to generate V m (D). For convenience, we embed it within a large circle, to eliminate unbounded edges. For ease of exposition, we will assume that every vertex has degree three10 . We will also assume that the cell containing each disc’s origin is known. Some terminology will be helpful. Let Vi− = V m (D) ∩ Pei . It has been shown ([41], Lemma 1) that this skeleton (e.g., the dotted lines of Figure 2.18) is connected. The core discs of Di are the set of discs that include Di and any discs intersecting Di . The core of Di is the union of cells of Di ’s core discs. For each edge e of V m (D), we determine if the possible cell for the site Di to its left needs constructing. If so, by Lemma 2.17, e is within Vi− , and every vertex of Pei can be found by performing, in a manner similar to that of ([35], Section 3.2), a depth first search (DFS) of Vi− starting from e. If the DFS is done in a systematic fashion, the vertices of Pei can be found in sorted order. As was the case in Section 2.2.3, a complication arises if the discs are not partially disjoint. If Di is a parent of Dj , then by Lemma 2.18, Rim will be empty; yet by Lemma 2.17, Pei has nonzero area. To deal with this problem, we define Dm to be the set of discs D where each disc’s radius ri has been replaced by rm − ri , rm being the radius of the largest disc of D. We will denote the modified version of each Di by Dim . It is easy to show that V m (D) is equivalent to the standard Voronoi diagram of the discs Dm , and that the discs of D that do not contain other 10 This assumption can be relaxed if a perturbation scheme (e.g., [35]) is employed. 41 discs (i.e., the ‘innermost’ discs) are exactly the outer discs of Dm (and the discs that contain these discs correspond to children within Dm ). Let us return to the situation where Di is a parent of Dj , and Dj itself does not contain any other discs of D. It follows that Djm is an outer disc of Dm , and Dim is one of its children. For each such child, we can (by Lemma 2.17) construct Pei by starting the DFS on any edge of the boundary of Rjm . As noted in Section 2.2.3, by using the algorithm of [26] to construct V m (D), we can determine the outer discs, and their corresponding sets of children, at no additional cost. Theorem 2.20 The possible Voronoi diagram of n uncertain discs can be conP structed in time O((n + k) log n + ni=1 g(i)), where k is the complexity of the diagram, and g(i) is the number of discs intersecting disc Di . Proof. We first use the procedure described above to construct V m (D) and to generate the possible regions. Next, we perform a plane sweep to construct the arrangement of these regions. We need only prove the running time. V m (D) can be constructed in O(n log n) time. For each Di , the DFS covers every edge in Vi− . We partition these edges into three sets: non-core edges, core boundary edges, and interior core edges. Note that non-core edges do not contain a cycle, else Pei would contain a non-core cell, contradicting Lemma 2.19. These edges thus form a set of 3-ary trees, each rooted at core boundary vertex, with vertices of Pei as leaves; hence the number of non-core edges is at most |Pei |, the complexity of Pei . For each core boundary vertex, there is a core boundary edge, so these number at most |Pei |. Since the interior core edges are a subset of a Voronoi diagram of the core discs, there are O(g(i)) of them. The total edges examined in the DFS is thus O(|Pei | + g(i)). P Since ni=1 |Pei | = O(k), generating all n possible cells can be done in time P O(n log n + k + ni=1 g(i)). The arrangement of the cells can then be generated, using a plane sweep, in time O(k log n), for the stated total running time. Theorem 2.21 The number of edges in a possible Voronoi diagram of n uncertain discs is O(n3 ), and this bound is tight. Proof. As noted earlier, Ve {} (D) is the arrangement of boundaries of possible cells Pei∈1...n . Its complexity is thus (i) the sum of the complexities of the individual 42 possible cells, plus (ii) the number of proper intersection points of edges from pairs of cells (Pei , Pej ). Since each possible cell is a cell of an additively weighted Voronoi diagram of n sites (which has O(n) complexity [66]), (i) is n·O(n). We can bound (ii) by noting that each point p that lies on the boundaries of Pei and Pej is the center of two discs: Ci , which has outside-tangent Di and inside-tangent some Dk , and Cj , which has outside-tangent Dj and inside-tangent Dm . Since no discs of D can lie in the interior of Ci or Cj , it follows that Ci = Cj . We can assume, without loss of generality, that rj > 0; for if ri = rj = 0, Pei cannot properly intersect Pej . This implies that j 6= k, since Dj cannot be both inside- and outside-tangent Ci . Each point p is therefore determined by an ordered triple of distinct discs (Di , Dj , Dk ). The total complexity of Ve {} (D) is thus O(n · O(n) + n ) = O(n3 ). 3 We can show that this bound is tight; see Figure 2.19. We place two sets of n/3 discs of large radii r clustered tightly around the points (−r − , 0) and (r + , 0), for small positive . We place an additional n/3 discs with zero radii in a stack near the origin. For each large disc Di , this stack induces Ω(n) arcs in the boundary of Pei ; and for each pair of large discs (Di , Dj ) on opposite sides of the x-axis, the boundaries of Pei and Pej will intersect in Ω(n) points. The total number of intersections in Ve {} (D) is thus (n/3)2 · Ω(n) = Ω(n3 ). Figure 2.19: Ve {} (D) complexity for discs is Ω(n3 ) (some edges omitted for clarity). 43 2.4 Guaranteed Delaunay Edges In this section, we consider the Delaunay triangulation of a set of sites, where each site lies within a region of uncertainty. Our goal will be to find an algorithm to determine the guaranteed Delaunay edges for these uncertain regions. Such an algorithm would be useful, since it would allow us to preprocess a set of these regions to determine a set of Delaunay edges that is a subgraph of every feasible Delaunay triangulation of the regions. Having this subgraph pre-computed could save a significant amount of the work required to (later) construct a Delaunay triangulation of a feasible set of the regions. Earlier, we defined a guaranteed Delaunay edge (a, b) to exist iff (a, b) is an edge in the Delaunay triangulation of every feasible set S = {s1 ∈ D1 , . . . , sn ∈ Dn } of D. This definition is not, however, as precise as we would like. If four sites of a set are cocircular, then their respective Voronoi cells will meet at a vertex of degree four. It is not immediately clear whether the cells of two sites that have only this vertex in common are considered to share an edge; what is more, the Delaunay triangulation of such a set is not unique (Figure 2.20). Figure 2.20: Cocircular points with two valid Delaunay triangulations (Voronoi diagram drawn with dotted lines). To deal with this problem, we will require guaranteed Delaunay edges to correspond to Voronoi edges of nonzero length. Formally, we will say that (a, b) is a guaranteed Delaunay edge of D iff for every feasible set S of D, the Voronoi cells for sa and sb share an edge of nonzero length. Our approach in determining the guaranteed Delaunay edges for a set of uncertain regions D will be to start with a set of pairs of these regions, then examine each 44 candidate pair (Da , Db ) to see if it corresponds to a guaranteed Delaunay edge. We S will let D denote the set Di ∈ (D \ {Da , Db }) . Lemma 2.22 If D intersects the convex hull of Da and Db , then (a, b) is not a guaranteed Delaunay edge. Proof. Suppose point q ∈ D ∩ CH(Da ∪ Db ) exists. First observe that as a consequence of the regions’ valid intersection property, there must exist some such q that does not lie on the boundary of either Da or Db . Now note that there exist points sa ∈ Da and sb ∈ Db such that q lies in the interior of sa sb (see Figure 2.21). This implies that there exists a feasible set of D in which q prevents sa and sb from CH(Da ∪ Db ) Da sa q sa q Db sb sb q sb sa Figure 2.21: Lemma 2.22. being neighbors in the Voronoi diagram of the set. Hence, (a, b) is not a guaranteed Delaunay edge. If Da and Db are regions of D, then we define a guaranteed Delaunay disc for (a, b) to be a disc that contains Da and Db , and does not intersect D. Lemma 2.23 If each region of D is convex, then (a, b) is a guaranteed Delaunay edge if and only if a guaranteed Delaunay disc exists for (a, b). Proof. If Cp is a guaranteed Delaunay disc for (a, b), then for every feasible set S of D, sa and sb lie within Cp , and Cp does not contain any other points of the set. We can shrink Cp until (without loss of generality) sa is on its boundary, then shrink it again, with respect to sa , to reach a disc Cp0 that also has sb on its 45 boundary. Hence, p0 lies on the Voronoi edge bordering cells of sa and sb ; and since Cp0 ⊆ Cp , no other sites of S intersect Cp0 , which implies that this Voronoi edge has nonzero length. As this holds for every feasible set S, (a, b) is by definition a guaranteed Delaunay edge. Now suppose (a, b) is a guaranteed Delaunay edge. For ease of exposition, we will assume that the uncertainty regions are discs, though the following proof can easily be adapted for any convex uncertainty regions. If Da and Db are not partially disjoint, then without loss of generality, Da contains Db . By Lemma 2.22, D does not intersect Da ; hence Da is a guaranteed Delaunay disc for (a, b). Suppose instead that Da and Db are partially disjoint. Consider the sequence of discs F that have inside-tangent both Da and Db , ordered by their origin’s position along hha, bii. We now prove that some element of F does not intersect D, and hence is a guaranteed Delaunay disc for (a, b). We start with an arbitrary disc Cp ∈ F . Let α (resp., β) be the tangent point of Cp with Da (resp., Db ), and let R (resp., L) be the portion of Cp that lies strictly to the right (resp., left) of directed segment αβ (Figure 2.22). Cp Da L α Db β R Figure 2.22: Lemma 2.23. If D intersects neither R nor L, then we are done (by Lemma 2.22, D cannot intersect αβ). Otherwise, assume by way of contradiction that D intersects both R 46 and L. Since every element of D is convex, the same Dk ∈ D cannot intersect both R and L (without intersecting αβ); hence, distinct discs Di and Dj of D intersect R and L respectively. It follows that there exists a feasible set S where sa = α, sb = β, si ∈ R, and sj ∈ L. We can apply Lemma 2.13 to this set to show that the Voronoi edge for (a, b) either consists of a single point (if D does not penetrate Cp ), or does not exist (otherwise). In either case, (a, b) cannot be a guaranteed Delaunay edge, a contradiction. We can then assume, without loss of generality, that only R intersects D. Let us now replace Cp with its neighbor within F whose origin lies to the left of αβ, and repeat the above procedure. This process must terminate, since at some point before Cp becomes infinitely large (and R = ∅), neither R nor L will intersect D. Interestingly, Lemma 2.23 does not hold for non-convex uncertainty regions. For example, for the regions shown in Figure 2.23, no guaranteed Delaunay disc exists for (a, b); yet for every feasible set of the regions, the Voronoi cells of Da and Db will share a nonzero-length edge, which implies that (a, b) is a guaranteed Delaunay edge. Dc Da Db Figure 2.23: Lemma 2.23 does not apply for nonconvex regions. From this point on, we will assume that D is a set of uncertain discs. We will make use of V 2 (D), the standard order-2 Voronoi diagram of the discs. Recall that 47 we denote the cell of this diagram corresponding to a disc pair (Di , Dj ) by R{i,j} . The spine of a pair of uncertain discs (Da , Db ) (denoted sp(a, b)) is the set of origins of the guaranteed Delaunay discs for (a, b) that have inside-tangent both Da and Db . Spines will be fundamental to our algorithm. We can identify some of their properties: Lemma 2.24 If (Da , Db ) is a partially disjoint pair of uncertain discs, then (a, b) is a guaranteed Delaunay edge if and only if sp(a, b) is nonempty. Proof. Assume (Da , Db ) is a partially disjoint pair as stated. If there exists a point p ∈ sp(a, b), then p is the center of a guaranteed Delunay disc for (a, b); hence, by Lemma 2.23, (a, b) is a guaranteed Delaunay edge. If (a, b) is a guaranteed Delaunay edge, then there exists some guaranteed Delaunay disc Cp for (a, b). We can shrink Cp until it has inside-tangent one of Da or Db , then shrink it further with respect to that tangent point until it has insidetangent both Da and Db . Since the final disc is within Cp , it must be a guaranteed Delaunay disc; hence, its origin lies on sp(a, b). Observe that if the spine of a pair of uncertain discs is nonempty, the pair must be partially disjoint. Lemma 2.25 If (Da , Db ) are partially disjoint discs from D, and sp(a, b) is nonempty, then R{a,b} is a nonempty cell within V 2 (D). Proof. Let p be a point on sp(a, b). The disc centered at p that has inside-tangent both Da and Db does not intersect D; hence, p is strictly closer to both Da and Db than to any point within D, which implies that p is in the interior of R{a,b} . Lemma 2.26 sp(a, b) is a connected subset of a hyperbolic arc. Proof. Every point p ∈ sp(a, b) satisfies d(p, oa ) + ra = d(p, ob ) + rb , so p lies on a hyperbolic arc. We now prove that the spine is connected. Place the axes so oa and ob are on the x-axis, with oa to the left of ob . Each point on the 48 spine now has a unique y-coordinate. Assume by way of contradiction that (i) u, v, and w are points on the hyperbolic arc containing sp(a, b), (ii) uy < vy < wy , and (iii) u, w ∈ sp(a, b), and v ∈ / sp(a, b). Let discs Cu ,Cv ,Cw be the discs centered at u, v, w that have inside-tangent Da and Db . Since v ∈ / sp(a, b), some point q ∈ Dk∈{a,b} ∈ D exists that penetrates / Cv . Let L be the line through the points of tangency of Cv with Da and Db . Observe that the two intersection points of the boundaries of Cu and Cv lie on or above L, while those of Cv and Cw lie on or below L; so if q lies on or above L, it penetrates Cw , and if it lies below L, it penetrates Cu . Thus either u or w is not in sp(a, b); a contradiction. Hence sp(a, b) is connected. We are ready to describe a procedure for constructing the guaranteed Delaunay edges of a set of n partially disjoint uncertain discs. We start by constructing V 2 (D), the order-2 Voronoi diagram of the discs. By Lemma 2.25, only those discs (Da , Db ) which have a nonempty cell R{a,b} within this diagram can be guaranteed Delaunay edges. Let us examine how the spine corresponding to each region can be efficiently constructed. If p is an endpoint (not at infinity) of sp(a, b), there must exist a disc Cp that has inside-tangent both Da and Db , and has outside-tangent some disc Dk contributing to D. By the proof of Lemma 2.25, p ∈ R{a,b} . We can shrink Cp with respect to Dk until it reaches disc Cq that intersects both Da and Db , and (without loss of generality) has outside-tangent Da . Note that q lies on the boundaries of both R{a,b} and R{a,k} . It follows that every such endpoint p of sp(a, b) is induced by a disc associated with a neighboring cell of R{a,b} , and we can construct the spine as follows. We initialize the spine to be the complete hyperbolic arc hha, bii. The we clip the spine to the hyperbolic arcs containing hha, kii (resp., hhb, kii) for each neighboring cell R{a,k} (resp., R{b,k} ) of R{a,b} . If the spine is nonempty after this process, then by Lemma 2.24, (a, b) is a guaranteed Delaunay edge. V 2 (D) contains O(n) edges, and can be generated in O(n log n) time [62]. The fact that the spine is connected (Lemma 2.26) implies that each clipping operation can be performed in constant time, and there are at most two such operations for each edge in V 2 (D); thus the running time is dominated by the O(n log n) time 49 spent constructing V 2 (D). This procedure can be modified to handle discs that may not be partially disjoint. As described in Section 2.2.3, we can determine (in O(n log n) time) which discs of D are outer discs, and which are the parents of those that are not. By Lemma 2.23, if Da is an outer disc with a child Db , then the only guaranteed Delaunay edge involving either Da or Db is (a, b), and only if no third disc intersects Da (note that this implies that Db is the only child of Da ). Lemma 2.27 If D is a set of partially disjoint discs, and Da , Db ∈ D intersect, then Da must be Voronoi neighbors with at least one disc of D that intersects Da . Proof. If Ra and Rb are neighboring cells in the Voronoi diagram of D, then we are done. Otherwise, let p be the intersection of oa ob with bisector hha, bii. We have d(p, oa ) − ra = d(p, ob ) − rb . Since Da and Db intersect, d(p, oa ) + d(p, ob ) ≤ ra + rb , and the two equations can be combined to show that d(p, oa ) ≤ ra , or in other words, p ∈ Da . Now, some point q in the interior of oa p must lie on bisector hha, cii, where Dc is a neighbor of Da . Hence d(q, oc ) − rc = d(q, oa ) − ra ≤ d(p, oa ) − ra ≤ ra − ra , which implies that q ∈ Dc . Our modified algorithm consists of these steps: 1. Construct V (D), the (order-1) Voronoi diagram of D. As a side effect, construct S, the subset of outer discs of D. Observe that since only discs of S will have cells within V (D), V (D) equals V (S). 50 2. Construct V 2 (S), the order-2 Voronoi diagram of S. 3. For each Da ∈ S with exactly one child Db , examine the neighboring discs of Da within V (D) to see if any other discs intersect Da . If none are found, then (by Lemma 2.27) Da is a guaranteed Delaunay disc for (i, j); output (a, b) as a guaranteed Delaunay edge. 4. For each cell R{a,b} of V 2 (S), if either Da or Db has children, do nothing. Otherwise, construct sp(a, b) as in the original algorithm, and output (a, b) as a guaranteed Delaunay edge if sp(a, b) 6= ∅. Steps 1 and 2 can be performed in O(n log n) time ([26],[62]). Steps 3 and 4 involve linear traversals of the cells of the two Voronoi diagrams. Each cell is traversed at most once, hence the steps take O(n) time. Since the total running time is the same as the original algorithm, we can claim: Theorem 2.28 The guaranteed Delaunay edges of n uncertain discs can be found in O(n log n) time. 2.5 Guaranteed Gabriel Edges of Uncertain Discs In this section, we present an algorithm to find the guaranteed Gabriel edges of a set of uncertain discs. As was the case for Delaunay triangulations, it is reasonable to suppose that constructing such a set as a preprocessing step could be useful for any application that required the Gabriel graph of a feasible set to be computed later. 2.5.1 Gabriel Discs and Zones From this point on, we will assume that D is a set of uncertain discs. We define a Gabriel disc of a pair of discs Da , Db ∈ D to be a disc with diameter αβ where α ∈ Da and β ∈ Db . We denote this Gabriel disc by Gα,β . The Gabriel zone of the pair is the union of all Gabriel discs of the pair, and is denoted Za,b . Formally, S Za,b = α∈Da ,β∈Db Gα,β . An example of a Gabriel zone is shown in Figure 2.24. 51 Za,b Db Da Figure 2.24: The Gabriel zone for a pair of uncertain discs. Lemma 2.29 (a, b) is a guaranteed Gabriel edge of D iff Za,b ∩ D = ∅. Proof. If Za,b ∩ D = ∅, then for every feasible set S of D, Gabriel disc Gsa ,sb will not intersect D; hence (a, b) is a guaranteed Gabriel edge. If there exists a point q ∈ Za,b ∩ D, then there exists some Gabriel disc Gα,β containing q. It follows that (a, b) is not a guaranteed Gabriel edge. Lemma 2.30 The Gabriel zone of discs Da and Db contains the convex hull of the discs. Proof. Every point in CH(Da ∪ Db ) lies on a segment αβ for some α ∈ Da , β ∈ Db ; and each such segment is the diameter of Gabriel disc Gα,β . Lemma 2.31 Every boundary point of Za,b is the source for a pair of perpendicular rays defining a quarter plane, where Da is tangent one ray, Db is tangent the other, and the disc centers lie within the quarter plane. Proof. Let p be a boundary point of Za,b . Since Za,b is a union of Gabriel discs, p must lie on the boundary of some Gabriel disc Gα,β where θ(∠αpβ) = 52 π 2. Con- → → and − sider the quarter plane bounded by the rays − pα pβ. Suppose, by way of contradiction, that p does not satisfy the stated properties (i.e., one of the discs is not tangent to one of the rays, or its center is outside the quarter plane). Without loss → and of generality, we can then assume that the quarter plane lies to the left of − pα, → (Figure 2.25). Now, consider the line L Da penetrates the right half-plane of − pα p0 p L α0 α β Da Db Figure 2.25: Lemma 2.31. → and tangent D at point α0 . Let p0 be the that is parallel to and to the right of − pα, a − → intersection of L with the line containing pβ. Observe that p0 lies on the boundary of Gabriel disc Gα0 ,β , and p is strictly interior to this disc. Hence p cannot lie on the boundary of Za,b . We will now investigate the structure of Gabriel zones in more detail. It will simplify our analysis if we place the coordinate axes so that Da is centered at the origin, and Db is centered on the positive x-axis. Let p be a point on the boundary of Za,b , and α, β be the points of tangency of the associated rays with the two discs. We will say that p is left-oriented (resp., → right-oriented) if β lies to the left of (resp., right of) − pα. Lemma 2.32 Segment oa ob lies in the kernel of Za,b . 53 Proof. Let K be the kernel of Za,b . We need only show that oa ∈ K, since a symmetric argument can be used to show that ob ∈ K, and the convexity of K then implies that oa ob ⊂ K. Let p be any point on the boundary of Za,b . By Lemma 2.31, p lies on the boundary of a Gabriel disc with diameter αβ, where pα is tangent to Da , and pβ is tangent to Db (Figure 2.26). Consider a point s moving from α toward oa . For each s, we construct disc C s , the circumcircle of points s, s0 , and β, where s0 is − → the projection of s onto ray pβ. Since ∠ss0 β is a right angle, sβ is a diameter of C s ; hence, C s is a Gabriel disc of Da and Db . Now observe that ss0 sweeps out a rectangle, one of whose diagonals is poa ; hence, poa ⊂ Za,b . As this holds for any boundary point p of Za,b , oa ∈ K. p s0 β α s ob oa Da Db Figure 2.26: Lemma 2.32. Lemma 2.33 Every boundary point of Za,b that lies to the left (resp., right) of oa ob is left-oriented (resp., right-oriented). Proof. We will first show that any points where the orientation of the boundary points change must lie on the x-axis (the line containing oa ob ). We will then show that this implies that there are exactly two such points. Finally, we will show that at least one point above (resp., below) the x-axis is left- (resp., right-) oriented. The proof then follows from the fact that the boundary is connected. 54 Let p be any transition point: a boundary point of Za,b that marks a transition → → and − from left- to right-orientation. Now, there must exist rays − pα pα0 tangent to −→ − → Da , and corresponding rays pβ and pβ 0 tangent to Db , as shown in Figure 2.27. Note that θ(∠αpβ 0 ) must equal θ(∠βpα0 ). This implies that p is collinear with the α0 Da β Db p ob β0 α Da Figure 2.27: Lemma 2.33: there are exactly two transition points. disc centers. Lemma 2.32 then implies that there are exactly two transition points, both lying on the x-axis. Now consider the vertical line that is tangent to Da at α and the horizontal line that is tangent to Db at β depicted in Figure 2.28. Observe that the intersection of these lines, p, lies within Za,b . Suppose the vertical line contains a right-oriented boundary point p0 of Za,b that lies above α. Clearly, its tangency point with Da −→ must equal α. By Lemma 2.31, the disc centers must lie to the right of p0 α; but this means p0 lies below α, a contradiction. Hence, every boundary point of Za,b on or above p (in fact, p is the only such point) must be left-oriented. A symmetric argument can be used to show that boundary points below the x-axis are right oriented. 55 β p Db Da α ob oa Figure 2.28: Lemma 2.33: boundary points above the x-axis are left-oriented. We will now derive an expression for points on the boundary of Za,b . Since Za,b is symmetric with respect to the x-axis, it will suffice to consider only upper boundary points: points above this axis. Let q be the common endpoint of perpendicular segments whose other endpoints are oa and ob respectively, and let θ = θ(∠ob oa q). Observe that q is one corner of a rectangle with sides of length ra and rb . Lemma 2.33 implies that p, the opposite corner of this rectangle, is a boundary point of Za,b ; see Figure 2.29. p ra q Da θ rb Db h Figure 2.29: Deriving an expression for boundary points of Za,b . 56 Let h be the distance between oa and ob . Then qx = h cos2 θ qy = h cos θ sin θ , and px = qx + (rb cos θ − ra sin θ) py = qy + (rb sin θ + ra cos θ) ; hence, px = h cos2 θ − ra sin θ + rb cos θ py = h cos θ sin θ + ra cos θ + rb sin θ . Gabriel zone boundaries are essentially generalizations of a class of curves known as pedal curves [44]. The pedal of a curve C with respect to a pedal point s is the set of points that are the intersections of two perpendicular lines, where one line is tangent to C, and the other contains s. A circle pedal curve results if C is a disc boundary. Observe that if the radius of Db is zero, then the (upper) boundary of the Gabriel zone for Da and Db is the circle pedal curve of the boundary of Da with pedal point s = Db . In the general case, rb > 0, and each point on the Gabriel zone boundary is generated by a line tangent to Da , and a second line that is tangent to Db instead of intersecting a fixed point. We call this generalization of pedal curves a bipedal curve, and when the curves are circles, as in our case, we call it a circle bipedal curve. To determine if some Dk contributing to D intersects the Gabriel zone Za,b , we proceed as follows. First, we transform Dk to the coordinate system of Figure 2.29. If Dk ’s origin is below the x-axis, then we replace it with its reflection through this axis (as noted earlier, Gabriel zones are symmetric about this axis, so Dk will intersect the zone if and only if its reflection does as well). Next, we find p, the nearest upper boundary point to ok , and calculate the directed line T tangent to the boundary at p (Figure 2.30). Finally, we return true if Dk intersects the left 57 half-plane of T . T ok p Dk Db Da Za,b Figure 2.30: Determining if Dk intersects Za,b . We will denote the coordinates of ok by kx and ky respectively. To find p, we will find the parameter θ that minimizes the square of the distance of boundary point p from ok , by setting the derivative of this distance to zero: d (px − kx )2 + (py − ky )2 = 0 dθ After algebraic manipulation, this reduces to finding the zeros of c1 cos θ sin θ + c2 sin2 θ + c3 cos θ + c4 sin θ + c5 58 where the coefficients are the constants c1 = h(2kx − h) c2 = 2hky c3 = ra kx − rb ky c4 = rb (kx − h) + ra ky c5 = −hky . The direction of the line T tangent to a boundary point p is: ∂px = −2h cos θ sin θ − ra cos θ − rb sin θ ∂θ ∂py = h(2 cos2 θ − 1) − ra sin θ + rb cos θ , ∂θ which simplifies to ∂px = −h cos θ sin θ − py ∂θ ∂py = −h sin2 θ + px . ∂θ Clearly, Dk will intersect the left half-plane of T iff the signed distance of ok from T is not more than rk . 2.5.2 Finding Guaranteed Gabriel Edges Our algorithm for generating the guaranteed Gabriel edges of a set of discs will consist of two steps. In the first, we generate a set of candidate edges; and in the second, we construct the Gabriel zone for each candidate and use it to determine whether the candidate is a valid guaranteed Gabriel edge. For the time being, we will assume that the uncertain discs are partially disjoint. Later, we will modify our algorithm to relax this restriction. Since the Gabriel graph of a set of points is a subset of the Delaunay triangulation of the points, the set of guaranteed Gabriel edges for D must be a subset of the guaranteed Delaunay edges of D. Hence, for the first step, we can use the same 59 set of candidates that we used for generating the guaranteed Delaunay edges in Section 2.4. Specifically, we can apply Lemma 2.24 and Lemma 2.25 to show that we need consider only candidates (a, b) that correspond to nonempty cells within V 2 (D). Let us now consider the second step of our algorithm. By Lemma 2.29, if (a, b) is not a guaranteed Gabriel edge, then some disc contributing to D must intersect Za,b . We will say that Dk is a certificate that (a, b) is not a guaranteed Gabriel edge. Determining whether (a, b) is a guaranteed Gabriel edge is thus equivalent to determining if no certificate for (a, b) exists. A simple approach for this problem would be to test each disc contributing to D for intersection with Za,b . This would require Ω(n) time for each candidate edge (a, b). In the following discussion, we show that we can improve this running time by reducing the number of uncertain discs that must be examined for each candidate edge. We will show that if a certificate for (a, b) exists, then some certificate exists that intersects a disc associated with a neighbor of Voronoi cell R{a,b} . Note that every cell adjacent to R{a,b} must be of the form R{a,x} or R{x,b} , where Dx is a disc contributing to D. We will denote this set of neighboring discs, {Dx }, by N . It is interesting to note that not every disc Dk ∈ D intersecting a Gabriel disc Cp is necessarily a member of N . In Figure 2.31, Dk is the only disc of D intersecting Gabriel disc Cp , yet N consists of the single element Dk0 (R{a,b} lies far below the discs shown in the figure). We will prove the following claim: Lemma 2.34 If D is a set of partially disjoint uncertain discs, cell R{a,b} is nonempty, and a certificate for (a, b) exists, then some member of N is a certificate for (a, b). Our proof will show that if no disc of N is a certificate for (a, b), then no certificate exists. We will do this by demonstrating that any Gabriel disc of Da and Db that intersects a certificate Dk can be manipulated to show that either Dk ∈ N , or that some surrogate certificate Dk0 ∈ N exists. We will require the following two lemmas. Lemma 2.35 If D is a set of partially disjoint uncertain discs, and R{a,b} 6= ∅, then R{a,b} contains a point on hha, bii. 60 Cp Dk β p α Da Db Dk0 Figure 2.31: Dk intersects Gabriel disc Cp , yet is not a member of N . Proof. Let p be a boundary point of R{a,b} (not at infinity). There exists a disc Cp and a disc Dk contributing to D where, without loss of generality (see Figure 2.32): (i) Cp is outside-tangent to Da and Dk , intersects Db , and is not penetrated by D; or (ii) Cp is inside-tangent both Da and Dk , is within Db , and is not interior to any disc contributing to D. If (i), then let ta be the point of tangency of Da with Cp . If ta ∈ / Db , then we can shrink Cp with respect to ta until it reaches a disc Cp0 that is outside-tangent to Db as well. Now, p0 ∈ hha, bii; and since D does not penetrate Cp0 , p0 ∈ R{a,b} . If instead, ta ∈ Db , then we can grow the zero-radius disc ta , maintaining its inside-tangency with Da , until it reaches a disc Cp0 inside-tangent Db as well. This must occur before it reaches Da , since the discs are partially disjoint. Observe that p0 ∈ hha, bii, and since ta is not in D’s interior, Cp0 is not interior any disc contributing to D. Hence, p0 ∈ R{a,b} . If (ii), we can grow Cp , maintaining its inside-tangency with Da , until it reaches disc Cp0 that is inside-tangent Db as well, so that p0 ∈ hha, bii. As with the previous case, Cp0 is not interior any disc contributing to D; hence, p0 ∈ R{a,b} . Lemma 2.36 If D is a set of disjoint uncertain discs, cell R{a,b} is nonempty, and Da ∩ Db ∩ D = 6 ∅, then some member of N is a certificate for (a, b). Proof. Suppose R{a,b} is nonempty, and I = Da ∩ Db ∩ D is nonempty as well. Observe that there must be a point p on the boundary of I that lies on the boundaries 61 Db Db Cp hha, bii Cp hha, bii p p Cp0 Cp0 0 ta Da Da Dk (i), ta ∈ / Db Da p ta Dk (i), ta ∈ Db Cp ta Dk Cp0 hha, bii Db (ii) Figure 2.32: Lemma 2.35. of (i) both Da and D, (ii) both Db and D, or (iii) both Da and Db . In the first two cases, p is a point on the boundary of R{a,b} , and intersects some disc Dk contributing to D; hence, Dk is within N , and (since p is a Gabriel disc of Da and Db ) Dk is a certificate for (a, b). In case (iii), p lies on hha, bii, at one of the two vertices of the lune Da ∩ Db . By Lemma 2.35, hha, bii ∩ R{a,b} 6= ∅; so starting with disc Cp = p, let us move Cp along hha, bii towards R{a,b} , maintaining its tangency with both Da and Db , until we reach Cq , where q is on the boundary of R{a,b} . If q ∈ I, then Cq is inside-tangent Da , Db , and some Dk contributing to D. Since q ∈ R{a,b} , Dk ∈ N ; and since Cq ⊂ I, Lemma 2.30 implies that Dk is a certificate for (a, b). If q ∈ / I, then Cq is outside-tangent Da , Db , and some Dk contributing to 62 D. Since q ∈ R{a,b} , Dk ∈ N . Note also that t, the point of tangency of Cq and Dk , must have appeared on the receding arc of Cp as it reached Cq , and that this receding arc lies within CH(Da ∪ Db ). Hence, by Lemma 2.30, Dk (which contains t) is a certificate for (a, b). We are now ready to prove Lemma 2.34. Suppose R{a,b} 6= ∅, and that Cp0 is a Gabriel disc (of Da and Db ) that intersects D. If Da , Db , and D have a common point of intersection, then we can apply Lemma 2.36, and we are done. Otherwise, since no common intersection point exists, we can shrink Cp0 until it first becomes outside-tangent to either Da , Db , or D. We can then shrink it further, with respect to this tangency point, until it reaches a disc Cp1 that intersects all three objects, and is outside-tangent two of them. Now, if one of the tangent objects is D, then p1 is on the boundary of R{a,b} , and the tangency point belongs to some Dk contributing to D, a member of N (Figure 2.33). In addition, since Cp1 ⊂ Cp0 , Dk is a certificate for (a, b). Dk Cp0 β p1 α p0 Cp1 Da Db Figure 2.33: p1 is on the boundary of R{a,b} . If the two tangent objects are Da and Db , then p1 lies on hha, bii (Figure 2.34). We can now use a similar approach to that of the proof of Lemma 2.36 to find some Cq that has outside-tangent Da , Db , and some Dk0 contributing to D. Clearly, Dk0 ∈ N . We need only show that Dk0 is a certificate for (a, b). There is one complication that did not occur in Lemma 2.36: t (the point of tangency of Dk0 with Cq ) still lies on the receding arc of Cq (shown in bold in Figure 2.34), but 63 this arc may not be contained by CH(Da ∪ Db ). Note however that the arc will lie within the union of CH(Da ∪ Db ) and the original Gabriel disc Cp0 . Hence, t lies within some Gabriel disc of Da and Db , and Dk0 (since it contains t) is a certificate for (a, b). Cp0 Cp1 Dk p1 α β p0 Dk0 Da Db t hha, bii Cq Figure 2.34: p1 lies on hha, bii. Now that we have proved Lemma 2.34, we are ready to describe an efficient algorithm for finding guaranteed Gabriel edges. Assume we are given D, a set of n partially disjoint uncertain discs. We first construct V 2 (D), in O(n log n) time, to generate the set of candidate edges. For each candidate (a, b), we examine each neighbor Dx of R{a,b} , to see if Dx intersects Gabriel zone Za,b . By Lemma 2.34, (a, b) is a guaranteed Gabriel edge iff no such intersecting neighboring disc is found. Since V 2 (D) has O(n) complexity, we will examine O(n) neighbors in total. Hence, the guaranteed Gabriel edges of D can be found in O(n · (C + log n)) time, where C is the time required to determine if a disc intersects a Gabriel zone. We can modify this algorithm to handle discs that are not partially disjoint. The approach is very similar to that employed in Section 2.4, and will make use of the following lemma. Lemma 2.37 Let S be the outer discs of uncertain discs D, and suppose D includes Da and Db , where Ra is nonempty, and Db is the only child of Da . If a certificate for (a, b) exists, then one of the neighbors of cell Ra ⊂ V (S) is a 64 certificate for (a, b). Proof. If Da intersects some Dk6=b ∈ D, then we can apply Lemmas 2.27 and 2.30 to show that some neighbor of Ra is a certificate for (a, b). Otherwise, if a certificate for (a, b) exists, there exists some Gabriel disc Cp0 (of Da and Db ) that intersects D. Since the diameter αβ defining Cp0 lies within Da , and Da ∩ D = ∅, we can shrink Cp0 until it is outside-tangent some Dk contributing to D (and not penetrated by D), at which time it will still intersect Da . We can then shrink the disc again, with respect to Dk , until it reaches some Cp1 that is outside-tangent Da and Dk . Now, Dk is a neighbor to Da in V (S); and since Cp1 ⊂ Cp0 , Dk is also a certificate for (a, b). Our modified algorithm consists of these steps: 1. We construct V (D), and S, the outer discs of D. 2. We construct V 2 (S). 3. For each Da ∈ S with exactly one child Db , we determine if any neighbors of Da within V (D) intersect Za,b . If none are found, then by Lemma 2.37, we output (a, b) as a guaranteed Gabriel edge. 4. For each cell R{a,b} of V 2 (S), if either Da or Db has children, we do nothing. Otherwise, as in the original algorithm, we examine the neighbors of R{a,b} to determine whether (a, b) is a guaranteed Gabriel edge. Since the modified algorithm has the same asymptotic running time as the original, we can claim: Theorem 2.38 The guaranteed Gabriel edges of n uncertain discs can be found in O(n · (C + log n)) time, where C is the time required to determine if a disc intersects a Gabriel zone. In Section 2.5.1, we showed that determining if a disc intersects Za,b can be performed efficiently using analytical methods (and a simple local coordinates transformation); hence, in practice, C = O(1), and the guaranteed Gabriel edges can be found in O(n log n) time. 65 2.6 Closing Remarks In this chapter, we provided complexity bounds and worst-case optimal algorithms for guaranteed Voronoi diagrams of uncertain discs and disjoint uncertain polygons. We introduced a new diagram, the possible Voronoi diagram of uncertain discs, and provided corresponding complexity bounds and an algorithm that is worst-case optimal (up to logarithmic terms) for its construction. We also demonstrated how to efficiently construct the guaranteed Delaunay and Gabriel edges of uncertain discs. One possible avenue for extending the results of this chapter would be to develop an algorithm for the guaranteed order-k Voronoi diagram of discs. Much of the theory we have presented concerning guaranteed (order-1) Voronoi diagrams can be applied to the order-k case, and an algorithm for constructing order-k diagrams might be very similar to that of Theorem 2.9. Perhaps the only complication with the order-k version is that no direct analog to Lemma 2.8 exists: the sequence of bisectors (hi, ji, . . .) forming the boundary of a guaranteed cell is not necessarily derived from a subset of the bisectors (hhi, jii, . . .) forming the boundary of the standard cell. To elaborate, let S be (the indices of) a subset of k discs, and p be a point on eS . As we did in for the order-1 case, we can the boundary of the guaranteed cell, R apply a mapping similar to that of Section 2.2.2 to find an analogous point q on the boundary of the standard cell, RS . Now, p will be the center of a disc Cp that has inside-tangent some Du∈S , and outside-tangent some Dv∈S / ; and q will be the center of a disc Cq that has outside-tangent some u0 ∈ S, and outside-tangent some 0 0 Dv0 ∈S / . As with the order-1 case, v will equal v; but unlike that case, u may be different than u. Basically, the disc of S that induces a point on the boundary of the guaranteed cell may be different than that which induces the point on the boundary of the standard cell. We can show that the inside-tangent discs u are closely related to the farthest disc Voronoi diagram of S, and hence can be efficiently determined by querying that diagram. What is more, we believe that an efficient algorithm for generating the order-k guaranteed Voronoi diagram can be found by extending the algorithm for order-1 diagrams and incorporating some of the techniques used in our 66 algorithm for constructing the guaranteed Voronoi diagram of uncertain polygons (Section 2.2.4), since it uses farthest point Voronoi diagrams to solve a very similar problem. Another possible extension to the results presented here would be to extend the results of Section 2.2.4 to allow the uncertain polygons to intersect. The (standard) Voronoi diagram of polygons is derived from that of line segments, and if the line segments are allowed to properly intersect, this has implications not only for the algorithm, but also for the complexity of the resulting Voronoi diagrams. Despite these concerns, we believe that the algorithm of Theorem 2.15 can be adapted with little or no changes to the case where the uncertain polygons are allowed to intersect, provided an efficient algorithm for constructing the standard Voronoi diagram of the polygons is available. Some other structures related to Delaunay triangulations that we did not consider, but that could be investigated within the framework of regions of uncertainty, are (Euclidean) minimum spanning trees and relative neighborhood graphs. In the relative neighborhood graph [71] of a point set, an edge exists between points a and b if no third point c exists that is closer to both a and b than they are to each other. These graphs were introduced by Toussaint [71], and he showed that they are subsets of Delaunay triangulations, and supersets of minimum spanning trees (in fact, they are subsets of Gabriel graphs as well). A natural extension of these graphs to the framework of uncertain regions would be to label two uncertain regions as guaranteed neighbors if they share an edge in the relative neighborhood graph of every feasible set of the regions. In its simplest formulation, the minimum spanning tree of a point set is the set of edges that (i) causes the points to be connected, and (ii) minimizes the sum of the lengths of the edges. One might define (a, b) to belong to the set of guaranteed spanning edges of a set of uncertain regions if (a, b) is an edge in the minimum spanning tree of every feasible set of the regions. We suspect that an efficient algorithm for generating such edges may be difficult to find, since unlike Delaunay, Gabriel, and relative neighborhood graphs, uncertainty involving sites in one part of the plane can affect whether two sites in a distant part of the plane are joined by a guaranteed edge. For example, in Figure 2.35, whether or not (a, b) is a guaranteed spanning edge depends upon the distance between c and d. 67 a c ? ? b d Figure 2.35: Guaranteed spanning edges. There is some existing research on minimum spanning trees of uncertain regions. Yang et al [75] showed that the problem of finding the feasible set of a set of uncertain discs or rectangles that minimizes the minimum spanning tree total edge length is NP-hard. Chambers et al [14] examined a related problem, that of finding the feasible set that minimizes the length of the minimum spanning tree’s longest edge. They showed that this problem is NP-hard for two classes of uncertain regions: unit length squares, and vertical line segments. A β-skeleton of a set of points is a generalization of a Gabriel graph, and was introduced by Kirkpatrick and Radke [39]. If p and q are two points in a set, then (p, q) is an edge in the β-skeleton iff for every third point r from the set, the angle ∠prq (or ∠qrp, whichever is smaller) is smaller than an angle derived from β. When β = 1, this angle is π 2, and the β-skeleton equals the Gabriel graph. For values β > 1, the analog to a Gabriel disc with diameter pq is a lune with vertices p and q. A natural extension of the results of this chapter to β-skeletons would start by stating that a guaranteed β-skeleton edge (p, q) exists iff for every feasible set S of a set of uncertain regions, the lune induced by p and q does not intersect S \ {p, q}. An applet exploring some of the geometric objects described in this chapter can be found at http://www.cs.ubc.ca/∼jpsember/gv.html. The applet (Figure 2.36) includes the following sample files: 1. discs.dat: Constructs the guaranteed Voronoi diagram of uncertain discs, 68 using an algorithm similar to the one described in the proof of Theorem 2.9. Additional objects that can be displayed include the possible Voronoi diagram of the discs, the possible cell of a particular disc, and the guaranteed Delaunay edges of the discs. 2. polygons.dat: Constructs a guaranteed Voronoi cell for one of a set of uncertain polygons. The program uses a simple, ‘brute-force’ algorithm, instead of that of Section 2.2.4. 3. gabriel.dat: Constructs the Gabriel zone for a pair of uncertain discs, and determines whether it is intersected by a third disc. Figure 2.36: Applet for Chapter 2. 69 Chapter 3 Convex Hull of Imprecise Points 3.1 Introduction Computing convex hulls is a fundamental problem in computational geometry. The concept of a convex hull is easy to grasp; imagine that the geometric plane is represented by a board, and a set of points in the plane is represented by a set of nails sticking out of the board. If a rubber band is stretched to enclose all of the nails and released, then the rubber band will come to rest in the shape of a convex polygon (Figure 3.1). A nail will be located at each vertex of the polygon, preventing the rubber band from contracting further. This polygon is the convex hull of the points. A set of points is convex if, for every pair of points a, b in the set, the segment ab is also within the set. The convex hull of a set of points S is the smallest convex set that contains S, and is denoted CH(S). Equivalent definitions for CH(S) include: the intersection of every convex set that contains S; and, the set of all points that can be expressed as convex combinations of points in S. There are many applications for convex hulls. In computer graphics, for instance, convex hulls can be used as an approximation of the boundary of an object. If no part of the convex hull of a shape is visible on the screen, it follows that no part of the shape is visible either. If the complexity of the shape’s convex hull is significantly less than that of the shape itself, then ensuring that the convex hull is at least partially visible before plotting the shape can lead to increased perfor70 Figure 3.1: The convex hull of a set of points mance. Convex hulls can also be used to calculate extrema associated with point sets. For example, consider the diameter of a set of points, which is the maximum distance separating any two points of the set. It is easy to show that the two points in question must lie on the convex hull of the points. Once the convex hull has been calculated, the technique of rotating calipers [72] can be applied to find these two points. Closely related to the diameter of a point set is its width, defined as the width of the smallest strip (or region bounded by parallel lines) that contains the set. This too can be found by applying rotating calipers to the set’s convex hull. Perhaps the simplest algorithm for constructing the convex hull of a set of points is the Jarvis march (also known as the gift wrapping algorithm) [33]. The algorithm constructs the ccw sequence of hull points, starting with the lowest point in the set (which is clearly guaranteed to lie on the hull). At each iteration, it calculates the next vertex by constructing a ray from the previous vertex through each of the original points. The point associated with the ray that makes the smallest ccw angle with the previous hull edge is output as the next hull vertex. It has an O(nh) running time, where n is the number of points, and h is the number of points on the convex hull. If every point lies on the convex hull (i.e., so that h = n), its running time is Θ(n2 ), which is suboptimal. O(n log n) algorithms for constructing convex hulls of points include Graham’s scan [29], and a divide-and-conquer algorithm from Preparata and Hong [58]. The Quickhull algorithm [9] is a randomized algorithm with expected O(n log n) 71 performance, and is so named due to its similarity to the Quicksort algorithm for sorting numbers. An output-sensitive algorithm is one whose running time depends not only on the size of the input, but also on the size of the output. An outputsensitive, O(n log h) convex hull algorithm was first proposed by Kirkpatrick and Seidel [40], and a simpler algorithm with the same running time was later given by Chan [15]. If the location of each point is unknown, but known to lie within a particular region of uncertainty, then the convex hull of the points (i.e., the ‘true’ convex hull) cannot be determined, and may be one of a continuum of feasible hulls: convex hulls of feasible sets. In this case, in order to construct a boundary approximation of the set, or to calculate extrema associated with the set, it may be useful to know which points in the plane are guaranteed to be within the true convex hull. To this end, we define the guaranteed hull of the set of uncertain regions as the set of points that are within the convex hull of every feasible set of the regions. If, instead, we wish to know which points in the plane may lie within the true convex hull, then we are interested in knowing the possible hull of the regions: the set of points that are within the convex hull of at least one feasible set. Figure 3.2 depicts a set of uncertain regions and its corresponding convex hull variants. The dotted boundary represents the convex hull of the regions, the dashed boundary represents the guaranteed hull of the regions, and the solid boundary represents the possible hull of the regions. As a motivating example for the possible hull problem, suppose each of the regions of Figure 3.2 is an island, and each island contains a sensor whose exact location is unknown. If a line-of-sight signal is transmitted between each pair of sensors, and the signal can detect any object passing between the pair, then a boat that wishes to approach the islands surreptitiously must remain outside of the possible hull of the islands. Possible hulls are also useful as a means of preprocessing a set of uncertain regions. Later, we will see how they can be used to identify those regions that contribute vertices to every feasible hull of the regions. 72 Figure 3.2: Convex hull, guaranteed hull, and possible hull of uncertain regions 3.1.1 Related Work An O(n log n) algorithm for constructing the guaranteed hull of a set of uncertain discs is given by Rappaport [61]. He shows that by modifying the discs’ radii in a particular way, the guaranteed hull of the original discs can be derived from the boundary segments of the modified discs’ convex hull. Using a domain-theoretic approach to computational geometry, Edalat et al ([21],[22]) define a partial convex hull as being a pair of open sets (E, I), where E is a portion of the plane guaranteed to be outside the hull, I is a portion guaranteed to be in the hull’s interior, and the hull boundary is exterior to both E and I. The closure of I is what we refer to as the guaranteed hull, while the closure of the complement of E is what we refer to as the possible hull. The authors prove that E and I can be calculated in O(n log n) time, when the uncertain input points are represented as partial points, which are essentially axis-aligned rectangles representing the uncertainty in the x and y coordinates. They extend this work to uncertain points represented by convex polygons in [23], and give an O(m · N log N ) algorithm for guaranteed hulls of such points, where N is the number of polygons, and 73 m is the number of distinct slopes of the polygons’ sides. Guaranteed and possible hulls of uncertain regions have also been investigated by Nagai, Seigo, and Tokura [55]. They provide O(n log n) algorithms for constructing the guaranteed hull of uncertain discs, and a particular class of polygon that includes axis-aligned rectangles. Nagai and Tokura expand this work in [54]; they represent a convex region as a function representing the distance of the object from a directed line with a particular polar angle, and show that the the convex hull of a set of regions represented in this way corresponds to the upper envelope of a set of such functions. Using this result, they achieve an O(n log n) algorithm for guaranteed hulls of uncertain convex polygons with n total vertices (an improvement over the algorithm of [23]). The same tools used to construct guaranteed hulls can be used to determine extremal feasible values for the width and diameter of imprecise points ([61],[55],[54]). For example, the minimum feasible diameter for a set of uncertain polygons with n total vertices can be found in O(n log n) time, while the maximum feasible diameter for the set can be found in expected O(nα(n) log n) time [54]. Löffler and van Kreveld have also studied convex hulls of uncertain regions [47]. They investigate the problem of optimizing a function among feasible hulls, such as finding convex hulls with maximal or minimal areas or boundary lengths, where the hulls are associated with regions of uncertainty that are line segments or axisaligned squares. They provide polynomial-time algorithms ranging from O(n log n) to O(n13 ) for some variants of these problems, and show that others are NP-hard. Variants involving uncertain discs are also discussed. Ezra and Mulzer [25] have investigated the problem of preprocessing uncertain regions to later allow efficient construction of convex hulls (when the uncertainty in the input has presumably been resolved). They show that with O(n2 ) preprocessing time, the convex hull of any feasible set of points drawn from a set of uncertain lines (or line segments) can be constructed in o(n log n) time. 3.1.2 Contributions This chapter looks at guaranteed and possible hulls of uncertain regions. We prove that the complexity of the guaranteed hull of a set of uncertain regions is strictly 74 linear in the number of the regions, regardless of the complexities of the individual regions. This improves the O(nα(n)) bound for uncertain polygons implicit in the paper by Nagai and Tokura [54]. We present two optimal algorithms for constructing guaranteed hulls of an important subclass of uncertain regions. This subclass includes discs and squares, the most common uncertain regions likely to be encountered in practice. The second of these algorithms, an adaptation of Chan’s convex hull algorithm [15], is outputsensitive (subject to some assumptions that we will give later). While Nagai and Tokura’s guaranteed hull algorithm [54] has a running time that matches that of our first algorithm, and can be modified so that its running time matches that of our second algorithm as well, it relies upon a dual-space approach that somewhat obscures the geometry of the uncertain regions. In our algorithms, this geometry remains explicit. In addition to these two new guaranteed hull algorithms, we show how three existing planar guaranteed hull algorithms can be modified to produce algorithms for 3-dimensional guaranteed hulls. We investigate possible hulls of nonconvex uncertain regions, and present an optimal algorithm for constructing such hulls for sets of uncertain polygons. We also show how possible hulls can be used within a preprocessing step to determine which regions contribute vertices to every feasible hull. 3.2 Guaranteed Hull The guaranteed hull of a set of uncertain regions D is defined formally as the set of points p such that for every feasible set S of D, p ∈ CH(S). We will denote the guaranteed hull of D by GH(D) (or, when D is clear from the context, simply by GH). One way in which guaranteed hulls can be useful is in preprocessing a set of imprecise points to optimize some later operation on (precise) feasible sets. For example, consider the problem of filtering out the regions (of a set of regions) that never contribute feasible hull vertices. If these regions can be identified, then their corresponding sites within any feasible set can be ignored when constructing feasible hulls. It is easy to show that a region has this property if and only if it lies 75 within the guaranteed hull of the other regions. The guaranteed hull of a set of uncertain regions is closely related to the set of directed lines that are tangent to pairs of the regions. This was shown in [55], though with a more complicated analysis that obscured this relationship to some extent. In this section, we will provide a geometrically intuitive definition of these tangent lines. For simplicity, we will assume that each set of uncertain regions D that we discuss has at least three regions. In addition, we assume that no line exists that is tangent to three regions of D, and also that no vertical line exists that is tangent to two regions of D. If a directed line L is right-tangent1 to a region of D, then we say that the region supports L. By contrast, if a region lies strictly to the right of L, then we say that the region invalidates L (Figure 3.3). A directed Di L Dj Figure 3.3: L is supported by Di , and invalidated by Dj . line is a hull tangent of D if some region of D supports the line, and the line is not invalidated by any other region of D. Observe that each angle in [0 2π) is the polar angle of a unique hull tangent of D. A bitangent of uncertain regions Di and Dj is a directed line L that is righttangent to both Di and Dj , with L ∩ Di preceding L ∩ Dj along L (by our general position assumption, these two intersections must each consist of a single point, and must be disjoint). Figure 3.4 shows some of the bitangents associated with some sets of uncertain regions. Note that our definition of bitangent differs from the standard definition, in which L may be left-tangent to either Di or Dj instead. We denote the bitangent of Di and Dj by Bij . For convenience, we will abbreviate θ(Bij ) (the polar angle of Bij ) as θij . No two bitangents induced by a set of uncertain regions can lie on the same line, as a consequence of our general position assumption. It is possible, however, 1 Recall that L is right-tangent to A iff (i) A intersects L, and (ii) A lies in L’s right half-plane. 76 Dj 0 Bkl Bkl Bij Di Dl Dk Figure 3.4: Bitangents for an ordered pair of regions to be associated with two (or more) bitangents (e.g., 0 of Figure 3.4). Bkl and Bkl We will say that Bij is a hull bitangent if it is also a hull tangent. By contrast, a hull tangent L is a pure hull tangent if it is not also a bitangent; e.g., exactly one region is right-tangent to L. There is a close connection between hull bitangents and guaranteed hulls: Theorem 3.1 If GH has nonzero area, then GH equals the intersection of the left half-planes of B, the hull bitangents of D; otherwise, either B = ∅, or the intersection of their left half-planes does not have finite, nonzero area. Proof. If X is a set of bitangents, we will denote the intersection of their left halfplanes by I(X ). We first consider the case where GH has nonzero area. start by proving that GH = I(B 0 ), where B0 We is the set of directed lines whose left half-planes intersect every region of D. If some point p ∈ GH is not within I(B 0 ), then p lies strictly to the right of some line L ∈ B 0 . Since every region of D intersects L’s left half-plane, there exists a feasible set S of D that lies to the left of L. This is a contradiction, since it implies that p ∈ / CH(S). Hence, p ∈ I(B 0 ), and GH ⊆ I(B 0 ). To show that I(B 0 ) ⊆ GH, let p be any point in I(B 0 ), and suppose p ∈ / GH. Then for some feasible set S, p ∈ / CH(S), which implies that S lies strictly left of some line Q through p. Let L be the hull tangent of D that is parallel to Q, and Di ∈ D be its right-tangent region. Since p ∈ I(B 0 ), L must lie on or 77 to the right of Q; but then S lies strictly left of L, a contradiction (since no point of Di lies to the left of P ). We now show that I(B 0 ) = I(B) by demonstrating that every boundary point s of I(B 0 ) lies on a line of B. Since I(B 0 ) and I(B) are convex, and I(B 0 ) ⊆ I(B), it then follows that I(B) = I(B 0 ) = GH. Let s be any boundary point of GH, and q be any interior point of GH. There must exist a line L0 ∈ B 0 through s; let L0 be the closest such line to q. If L0 is right-tangent to two regions of D, then L0 is a hull bitangent, and we are done. Otherwise, we can reach a contradiction as follows. If L0 is not right-tangent to any regions of D, we can translate L0 to the left to yield another line within B 0 that has s strictly to its right, a contradiction since this implies that s ∈ / GH. If L0 is right-tangent to a single region of D, then we can reach a contradiction by rotating L0 around the single tangency point to yield an element of B 0 that either has s strictly to its right, or (if s is the tangency point) is closer to q. Now let us consider the case where GH has zero area (which includes the case where GH = ∅). Suppose, for the sake of contradiction, that B = 6 ∅, and that I(B) has bounded, nonzero area. There must exist a point s that is interior to I(B) and exterior to GH. In addition, s must lie strictly to the right of L, some hull tangent of D. Let Di be the region right-tangent to L, t be its tangency point, and N be the line through s that is perpendicular to L (Figure 3.5). If we now now rotate L cw (or ccw, if t lies to N ’s right) around Di , it must either reach a bitangent of D, or become parallel to and strictly left of some bitangent containing an edge of B, before it reaches s. We now have a contradiction: the former implies that s ∈ / I(B), and the latter implies that the original L was not a hull tangent. Observe that since the set of lines right-tangent to a region is equal to the set of lines right-tangent to a region’s convex hull, Theorem 3.1 implies that each uncertain region can be replaced with its convex hull without affecting the guaranteed hull of the regions. When ordered by their polar angles, the hull tangents for a set of uncertain regions have a particular structure: Lemma 3.2 For every hull tangent L of D supported by a region Di , there exist adjacent bitangents (. . . B∗i , Bi∗ . . .) in the polar angle-ordered sequence of hull bitangents of D where θ(L) ∈ [θ∗i θi∗ ]. 78 N I(B) s Di L t Figure 3.5: Theorem 3.1. Proof. We prove the claim by showing that if L is a pure hull tangent, and is rotated ccw (resp., cw) around Di , the first hull bitangent it reaches will be of the form Bi∗ (resp., B∗i ). We will prove only the ccw case, as the cw case is very similar. First, observe that Di cannot contain any other region of D, since that region would invalidate L. Note also that since L is a pure hull tangent, every other region of D penetrates its left half-plane. Let M be the first hull bitangent reached when L is rotated ccw around Di , and pi ∈ Di and pj ∈ Dj be its associated tangent points. Suppose pj precedes pi on M . We can now rotate M slightly cw around Di to reach line N that is right-tangent to Di and has Dj strictly to its right (Figure 3.6). Now, Dj penetrated the left half-plane of (the original) L; hence Dj Di pj L N pi M Figure 3.6: Lemma 3.2. L must have become right-tangent to Dj before M . This is a contradiction, since 79 we would have stopped rotating L at that time. Thus, pj follows pi along M , and M = Bij . Note that a similar argument will show that if L is invalidated by some other region of D (and hence is not a hull tangent), then if it is rotated ccw (resp., cw) around Di , the first hull bitangent it reaches will be B∗i (resp., Bi∗ ). Besides their role in constructing guaranteed hulls, we can identify another use for hull bitangents: Theorem 3.3 If L is a directed line, B is the angle-ordered set of hull bitangents of uncertain regions D, and c is the time required to construct tangent lines for regions of D, then it can be determined in O(log |B| + c) time whether the left half-plane of L contains a feasible set of D. Proof. By Lemma 3.2, B will contain adjacent bitangents (. . . , Bij , Bjk , . . .) where Dj supports a hull tangent M that is parallel to L, and θ(L) ∈ [θij θjk ]. Note that the left half-plane of L will contain a feasible set of D iff M lies in the left half-plane of L. We can perform a binary search, in O(log |B|) time, to identify Dj ; and we can then construct M , in time c. We define an aligned similar set of uncertain regions to be a set with the property that each region in the set can be derived from a common base shape by applying only translation and scaling transformations. We say that each such region is an instance of this base shape. We further require such sets to have the following properties: 1. Each region is compact (i.e., closed and bounded). 2. It can be determined in constant time which of two instances of the base shape is larger. 3. In constant time, the line with a particular polar angle that is right-tangent to an instance of the base shape can be constructed. 4. In constant time, the bitangent Bij of instances Di and Dj can be constructed, or it can be determined that no such bitangent exists. 80 Aligned similar uncertain regions include sets of discs with varying radii and origins, and sets of axis-aligned squares. We will say that two points are similar if they belong to distinct instances of a common base shape, and correspond to the same point in the base shape. For example, the lower left corners of a pair of axis-aligned squares are similar points. We conclude this section by proving two properties of aligned similar uncertain regions that we will make use of later: Lemma 3.4 If Di and Dj are aligned similar uncertain regions of D, and Di and Dj are partially disjoint, then Bij and Bji are unique and distinct, and the polar angle of every hull tangent supported by Di is within [θji θij ]. Proof. Let Di and Dj be as stated, H be their convex hull, and S be the set of directed lines right-tangent to H (ordered by polar angle). Observe that we can partition S into a sequence of four subsets, where each subset contains, in order: lines right-tangent only to Di , a bitangent of the form Bji , lines right-tangent only to Dj , and a bitangent of the form Bij . It follows that Bij and Bji are unique and distinct, and that every line right-tangent to Di with polar angle within (θij is a member of the first subset (and is thus invalidated by Dj ). θji ) Lemma 3.5 If Di and Dj are aligned similar uncertain regions, Di and Dj are partially disjoint, and Di is no smaller than Dj , then [θji θij ] spans at most π radians. Proof. Let Di and Dj be two such regions. Place the axes so that Bij is aligned with the negative x-axis. Let pi ∈ Di , pj ∈ Dj be the tangent points of Bij , and qi ∈ Di , qj ∈ Dj be the tangent points of Bji ; see Figure 3.7. Since the function that transforms points of Di to the similar points of Dj preserves angles, pi and qi → −−→ are similar, as are pj and qj . This implies that θ(− p− i qi ) = θ(pj qj ). Note also that since Di is no smaller than Dj , d(pi , qi ) ≥ d(pj , qj ) > 0. Hence the angle that Bji makes with Bij (the negative x-axis) is in the range (0 π]. 3.3 Guaranteed Hull Complexity In this section, we show that the guaranteed hull of a set of n uncertain regions has O(n) complexity, regardless of the complexities of the regions themselves. 81 Bji qi qj Di Dj pi pj Bij Figure 3.7: Lemma 3.5 This improves the O(nα(n)) complexity bound for guaranteed hulls of uncertain polygons implied in [54]. Since hull bitangents are closely tied to guaranteed hulls, we also provide tight bounds for the number of hull bitangents induced by sets of aligned similar uncertain regions. We also discuss an open problem, the number of hull bitangents induced by uncertain polygons. Theorem 3.6 The guaranteed hull of n uncertain regions has O(n) complexity, and for some sets of uncertain regions, this bound is tight. Proof. Let D be a set of n uncertain regions. By Theorem 3.1, GH is the intersection of the left half-planes of the hull bitangents of D. We will show that at most n bitangents contribute to GH by showing that for a particular region Di ∈ D, at most one bitangent of the form Bi∗ contributes to GH. Suppose, by way of contradiction, that bitangents Bij and Bik contribute to GH. Without loss of generality, we can assume that θij < θik (Figure 3.8). Let pi ∈ Di be the tangent point of Di associated with Bij . Observe that if we rotate Bij slightly cw around pi , every region of D will intersect the left half-plane of the (rotated) bitangent. Hence, some feasible hull of D lies in this half-plane, which implies that no point on the original bitangent before pi lies within GH. Since Bij contributes to GH, there must exist a point x ∈ Bij that is a non-vertex boundary point of GH. Now, to be within GH, x must lie to the left of Bik ; but then x lies before pi , a contradiction. For some sets D, the bound is tight; for example, if D consists of n small discs each centered at the vertices of a large regular n-gon, GH(D) is a (slightly smaller) regular n-gon. 82 Dj Dk Di pi Bik Bij Figure 3.8: Theorem 3.6 Observe that the complexity of a guaranteed hull is linear in the number of uncertain regions, and is independent of the complexity of the individual regions. This improves the O(nα(n)) bound for guaranteed hulls of uncertain polygons given in [54] (p. 259, Corollary 11). While the complexity of the guaranteed hull is linear at worst, no linear bound has been proven for the number of hull bitangents induced by a set of uncertain regions. For discs at least, the number of hull bitangents is strictly linear ([55]). We can generalize this result: Theorem 3.7 If D is a set of n uncertain regions, and each pair of regions induces at most two bitangents, then (i) D induces at most 2n − 2 hull bitangents; and (ii) there exist sets D for which this bound is tight. Proof. Let D be one such set of n regions. The hull bitangents of D, ordered by polar angle, have the form (Bi1 ∗ , Bi2 ∗ , . . . , Bik ∗ ), which can be expressed compactly as the (cyclic) sequence (i1 , i2 , . . . , ik , i1 ). We will prove the upper bound by showing that this is a Davenport-Schinzel cyclic sequence [67] of order 2. Suppose, by way of contradiction, that (i) Da ∈ D induces hull bitangents Bau and Bav , (ii) Db ∈ D induces hull bitangents Bbw and Bbx , and (iii) θau < θbw < θav < θbx . Observe that by Lemma 3.4, θau and θav are both within (θba θab ], whereas θbw and θbx are outside of this range. Since this contradicts the ccw ordering of the bitangents, no such subsequence of hull bitangents exists. 83 The ordered hull bitangents thus form a Davenport-Schinzel cyclic sequence of order 2, the maximum length of which is 2n − 2 ([67], Corollary 1.10). This bound is tight, as Figure 3.9 illustrates.2 Figure 3.9: Lemma 3.7: 2n − 2 hull bitangents The bound of Theorem 3.7 applies to aligned similar uncertain regions, and also to any set of uncertain regions that are pairwise disjoint (since each pair of such regions will induce exactly two bitangents). If every region in a set of aligned similar uncertain regions is the same size, then we can tighten the bound: Lemma 3.8 If D is a set of n aligned similar uncertain regions of uniform size, then D induces at most n hull bitangents. Proof. Let S = {si ∈ Di } be any set of pairwise-similar points drawn from regions of D, and let L be the directed line containing si sj , where si and sj are any two points of S. Suppose we move L to its left until either (i) it reaches Bij , or (ii) some region of D lies strictly to its right. Since the regions are all translates of a common base shape, (i) will occur first if and only if si sj is an edge of CH(S). It follows that the hull bitangents of D correspond exactly to the edges of CH(S). The proof of Lemma 3.8 suggests the following approach for constructing guaranteed hulls of aligned similar uncertain regions of uniform size. We first select a set of similar points from the regions, e.g., by choosing each region’s leftmost point (assuming it is unique). We then construct the convex hull of these points. Any convex hull algorithm will suffice, as long as it produces the h hull vertices in ccw order and identifies which region contributed each vertex. We then construct the 2 Interestingly, while the regions of Figure 3.9 realize the upper bound, their guaranteed hull is actually empty. 84 h hull bitangents corresponding to these vertices, and calculate the intersection of their left half-planes (in O(h) additional time). The complexity of the guaranteed hull of a set of uncertain regions can be substantially smaller than the number of bitangents induced by the regions. In fact, it is easy to construct sets of uncertain regions that induce a linear number of hull bitangents yet have an empty guaranteed hull (e.g., Figure 3.9). For uncertain polygons with total complexity n, it was shown in [55] (using techniques from [6]) that the number of hull bitangents is O(nα(n)). It is not known whether this bound is tight. The polygons in Figure 3.10 suggest that it is unlikely that an argument using Davenport-Schinzel theory similar to that employed in the proof of Theorem 3.7 will work to improve this bound. In the figure, vertex a supports hull bitangents Bac and Bad , vertex b supports hull bitangents Bbe and Bbf , and θac < θbe < θad < θbf . Bac Dd Da Bbe Bad a De Db b Dc Bbf Df Figure 3.10: Some hull bitangents of a set of uncertain triangles 85 3.4 Guaranteed Hull of Aligned Similar Uncertain Regions In this section, we present two optimal algorithms for constructing guaranteed hulls of aligned similar uncertain regions. The second of these, an adaptation of Chan’s convex hull algorithm [15], is output-sensitive, if we consider its output to be the hull bitangents of the regions (and not the guaranteed hull, which is the intersection of their left half-planes). The first algorithm has the advantage of being simpler (and often faster in practice) than the second. The geometry of the uncertain regions remains explicit in both algorithms. This is in contrast to the dual-space approach used in Nagai and Tokura’s algorithm [54], which can be adapted to match the running time of our second algorithm as well. It will be convenient to focus our attention on the upper guaranteed hull of D, which is the intersection of the left half-planes of the hull bitangents whose polar angles lie in the range [ π2 , 3π 2 ]. Once the upper and (by symmetric procedures) lower hulls are found, the complete guaranteed hull is easily constructed as their intersection. For simplicity, throughout this section, we will assume that uncertain regions are disc-shaped (the algorithm can be easily adapted for general aligned similar uncertain regions). We define the east disc, De , to be the disc whose leftmost point is farthest to the right, and the west disc, Dw , to be the disc whose rightmost point is farthest to the left. Let B•e be the line right-tangent to De with polar angle π2 , and Bw• be the line right-tangent to Dw with polar angle 3π 2 . Observe that both B•e and Bw• are hull tangents of D. Note also that if Bw• is not strictly left of B•e , then GH is empty. We will assume from this point on that this is not the case (which implies that De 6= Dw ). 3.4.1 Algorithm 1 Our first algorithm has an O(n log n) running time. Starting with the east and west discs of a set of aligned similar uncertain discs, it maintains a polar angle-ordered sequence of upper hull bitangents for the set, and updates these bitangents as new discs are added to the set. This approach for generating the guaranteed hull is similar to the one used by 86 the standard convex hull algorithms of Shamos [64] and Preparata[57]. Those two algorithms are on-line algorithms: they process their input in a serial fashion, one item (in this case, a point) at a time, and must be able to maintain the solution for the inputs seen so far without delaying any processing until all inputs have been seen. The algorithm for uncertain discs that we will present is not on-line, since it will require a step that preprocesses the entire input; but once this preprocessing is complete, the algorithm behaves essentially as an on-line algorithm, since it maintains the list of upper hull bitangents corresponding to the set of discs for the subsequence of discs seen up to the current point. One surprising characteristic of this approach is that if the remaining discs (other than De and Dw ) are processed in an arbitrary order, a suboptimal algorithm must result. Consider the discs {A1 , B1 , B2 , B3 , De , Dw } of Figure 3.11. A1 B1 B2 B3 De Dw Figure 3.11: Hull bitangents of discs 87 The (upper) hull bitangents for these six discs are shown as directed lines. Each Bi induces a pair of hull bitangents with A1 . In fact, we could construct an example with discs {B1 , . . . , B(n−2)/2 }, inducing Θ(n) hull bitangents (the first two discs are De and Dw ). Now, consider a disc A2 , the same size as A1 , but whose origin is slightly above that of A1 , so that A2 lies strictly to the right of the existing hull bitangents. A2 would thus invalidate these bitangents, and induce a new set of (slightly different) hull bitangents. We can construct Θ(n) such discs, {A1 , . . . , A(n−2)/2 }, where each new disc Ai>1 invalidates Θ(n) hull bitangents of the previous set. Hence, any algorithm that incrementally maintains a list of hull bitangents as new discs are added, where the discs are processed in an arbitrary order, will require Ω(n2 ) steps in the worst case. In the remainder of this section, we will show how sorting the discs by size leads to an optimal algorithm. We will assume that D is the ordered set {De , Dw , D3 , . . . , Dn }, where the discs {D≥3 } are sorted by increasing size. We will denote by B the ordered sequence of upper hull bitangents maintained by the algorithm. To simplify our algorithm, we will include the vertical hull tangents B•e and Bw• in this sequence as well. At the heart of our algorithm is a binary search performed for each disc in the sorted sequence. We will show that if the disc contributes new hull bitangents (i.e., so that B changes as a result of the new disc), then a binary search can find the location where the new bitangents are to be added. An interesting point about the search is that if the new disc does not contribute new hull bitangents, then this fact can be deduced from the binary search’s output, whatever it ends up being. 1. Let D = {De , Dw , . . . , Dn }, where {D≥3 } are sorted by increasing size. 2. Initialize q to 2, and B to (B•e , Bew , Bw• ). 3. Increment q; if q > n, stop. 4. Find Bij , the last bitangent in B for which θiq ∈ (θij 3π 2 ). 5. If no such Bij is found, or Bij is not followed by Bjk where θjq ∈ (θij θjk ), then go to step 3. 88 6. Find Bmn , the first bitangent of B following Bij that is not invalidated by Bq . 7. Replace bitangents of B strictly between Bij and Bmn with the pair (Bjq , Bqm ), and go to step 3. Let us prove the correctness of the algorithm. We will denote the contents of B at the start of each step 3 by Bq . Lemma 3.9 If aligned similar uncertain discs Da and Db support hull tangents A and B respectively, where θ(A) < θ(B), then Bab exists and has a polar angle within the range [θ(A) θ(B)]. Proof. If Da and Db support hull tangents A and B, then they must be partially disjoint, which implies that Bab exists. Now, Lemma 3.4 implies (i) θ(A) ∈ [θba θab ], and (ii) θ(B) ∈ [θab θba ]. If θab < θ(A), then (i) implies that θ(A) ≥ θba > θab , violating (ii); and if θab > θ(B), then (ii) implies that θab > θba ≥ θ(B), violating (i). Hence, θ(B) ≥ θab ≥ θ(A). By Lemma 3.9, π 2 ); 3π 2 = θ(Bw• ) > θew > θ(B•e ) = π 2, while θwe ∈ ( 3π 2 hence Bew is the only upper hull bitangent of B2 , and step 2 is correct. The remaining steps of the algorithm are repeated for each additional disc Dq . The following lemma is a key to the algorithm’s efficiency: Lemma 3.10 Let q ∈ {3, . . . , |D|}. Then Bq is either Bq−1 , or can be constructed by replacing a contiguous subsequence of Bq−1 with a pair of bitangents (Bxq , Bqy ). Proof. If Dq contributes no bitangents to Bq , then Bq ⊆ Bq−1 ; and since Dq supports no hull tangents of {De , Dw , . . . , Dq }, no bitangent of Bq−1 is invalidated by Dq . Hence, Bq = Bq−1 , and we are done. Let us proceed assuming that Dq contributes bitangents to Bq . Suppose, by way of contradiction, that it contributes more than two. Then, by Lemma 3.2, it must contribute (at least) distinct hull bitangents of the form Bqi and Bqj , where θqi < θqj . Consider the case where Di is no larger than Dq . By 89 Lemma 3.4, θqj ∈ [θiq θqi ]; and by Lemma 3.5, θqi − π is a lower bound for θiq , and θiq +π is an upper bound for θqi . Hence, θqj ∈ [θqi −π θqi ]∩[θiq θiq +π], which (since Bqi is an upper hull bitangent) implies that θqj ≤ θqi , a contradiction. If instead, Di is larger than Dq , then Di is either the east or west disc. If Di is the east disc, then it supports B•e ; but since Di is larger than Dq , the above argument will show that π 2 ∈ [θqi θqi + π], a contradiction. If Di is the west disc, then Lemma 3.4 implies that θqj ∈ (θiq 3π 2 > θqj > θiq > θqi > that 3π 2 ∈ (θqi π 2. θqi ). Since θqj > θqi , this implies But since Di is the west disc, the same lemma implies θiq ), a contradiction. Thus, if Dq contributes any hull bitangents to Bq , it contributes exactly two: Bxq and Bqy . By Lemma 3.2, every hull tangent of Dq with polar angle θqy > θ > θxq must be a pure hull tangent supported by Dq . Hence, every bitangent of Bq−1 with polar angle in this range will not appear in Bq . All that remains is for us to show that every bitangent of Bq−1 with a polar angle outside of this range appears in Bq as well. Suppose Buv appears in Bq−1 but not in Bq , and that θuv > θqy . Now, Buv must have been invalidated by Dq , which implies that there exists a directed line L supported by Dq that is parallel to and strictly to the right of Buv . But then L is a hull tangent of Dq , so by Lemma 3.2, there exists in Bq some hull bitangent Bq∗ where θq∗ > θqy , a contradiction. A symmetric argument handles the case where θuv < θxq . Suppose we are starting step 4 of a particular iteration of the algorithm, so that Dq is the next disc to be processed. For the moment, assume Dq contributes upper hull bitangents to Bq . By Lemma 3.10, it must induce a new pair of hull bitangents (Bxq , Bqy ), and there may exist a subsequence of bitangents of Bq−1 that are invalidated by Dq . Bq−1 is the concatenation of three contiguous subsequences: 1. bitangents not invalidated by Dq whose polar angles are less than θxq ; 2. bitangents invalidated by Dq ; and 3. bitangents not invalidated by Dq whose polar angles are greater than θqy . Note that the first and third subsequences are nonempty, since they must include B•e and Bw• respectively. This implies that there must exist a last bitangent of the first subsequence. Let us now prove that this bitangent must be Bij of step 4. 90 Suppose Dq≥3 contributes upper hull bitangents to Bq , and that Bij is a bitangent in Bq−1 . We need to show that Bij belongs to the first subsequence iff 3π 2 ). θiq ∈ (θij Observe that Dq supports some pure upper hull tangent L of Dq . If Bij belongs to the first subsequence, then by Lemma 3.9, θ(L) > θiq > θij , 3π 2 ), which implies θiq ∈ (θij and we are done. If Bij does not belong to the 3π 2 ) first subsequence, then we can assume θiq ∈ (θij Hence, 3π 2 > θiq > θij > (otherwise we are done). π 2. If Bij belongs to the second subsequence, then since Dq invalidates Bij , Dq must support a hull tangent that is parallel to, and to the right of, Bij . It follows (by Lemma 3.4) that θij ∈ (θiq θqi ), which in turn implies that 3π 2 > θiq > θqi > π2 . Now, if Dq is at least as large as Di , then Lemma 3.5 implies that the range [θiq θqi ] spans at most π radians, a contradiction. If instead, Dq is smaller than Di , then Di must be either the east or west disc and support a hull tangent of Dq with polar angle π 2 or 3π 2 . Lemma 3.4 then implies that either π 2 or 3π 2 is within [θqi θiq ], a contradiction. Finally, if Bij belongs to the third subsequence, then since Bij is a hull bitangent of Dq , Lemma 3.4 implies that θij ∈ (θqi θiq ). Thus 3π 2 > θiq > θqi > π2 , and we can continue as in the previous case to reach a contradiction. Observe that Dq contributes upper hull bitangents to Bq if and only if consecutive bitangents (Bij , Bjk ) exist in Bq−1 , where θjq ∈ (θij θjk ) (i and k need not be distinct). As we have shown, if such a pair exists, then a binary search will return Bij in agreement with step 5. Now consider the situation where Dq does not contribute upper hull bitangents to Bq . As we mentioned earlier, regardless of what the binary search of step 4 returns, we can recognize that Bq remains unchanged. Since in this situation, no pair of consecutive bitangents (Bij , Bjk ) exists in Bq−1 that satisfy θjq ∈ (θij θjk ), one of three failure events must occur: the search will (i) fail to return a bitangent, (ii) return a bitangent with no successor (i.e., the last bitangent in the list), or (iii) return a bitangent Bij with successor Bjk where θjq 6∈ (θij θjk ). Since these are exactly the tests performed in step 5, if a failure event occurs, the algorithm continues with Bq = Bq−1 . Otherwise, steps 6 and 7 convert the bitangent sequence, from Bq−1 to Bq , in accordance with Lemma 3.10. Let us now determine the algorithm’s running time. Step 1 requires O(n log n) 91 time to determine the east and west discs, and sort the remaining discs by size. Step 2 can be done in constant time. To allow O(log n) insertions, deletions, searches, and traversals, we store the hull bitangents B in a balanced binary tree. As a consequence, each iteration of the remaining steps of the algorithm can be performed in O((1 + k) log q) time, where k is the number of bitangents removed3 from Bq . Since a bitangent can be removed only once, and at most two bitangents are added in each of the O(n) iterations, the total k over all iterations is O(n). Hence, the running time of the algorithm is O(n log n). One complication arises if the discs are not partially disjoint (i.e., some discs are contained by others). The key observation to make in this setting is that no disc that contains another can contribute to the discs’ guaranteed hull: Lemma 3.11 If Di , Dj ∈ D, and CH(Di ) contains CH(Dj ), then Di cannot contribute to GH(D). Proof. If CH(Di ) contains CH(Dj ), then by our general position assumption, every line right-tangent to Di is invalidated by Dj , and hence cannot be a hull tangent. The rest of the proof follows from Theorem 3.1. Suppose that in some iteration of the algorithm, Dq contains some other disc Dk of D. Since the discs are processed in order of increasing size, and the east and west discs never contain other discs, k must be less than q; it follows that Dq cannot contribute upper hull bitangents to Bq . When processing Dq , we thus expect our algorithm to encounter one of the three failure events described earlier. In fact, since no bitangent Bkq or Bqk exists, the algorithm may encounter a fourth type of failure event: an attempt to construct a bitangent that doesn’t exist. This additional test can be included in step 5. Theorem 3.12 The above algorithm can be used to generate the guaranteed hull of n aligned similar uncertain discs in O(n log n) time. Proof. As we have shown, the upper and lower guaranteed hulls can be found in O(n log n) time (note that the sorted disc sequence (De , Dw , . . . , Dn ) need be 3 If a suitable tree variant is used (e.g., [73]), each step can be done in O(log q + k) time. This will not change the asymptotic running time of the algorithm, however. 92 constructed only once). By Theorem 3.1, GH is equal to the intersection of the left half-planes of the upper and lower bitangents. Finding such an intersection is a basic problem of computational geometry, and an optimal O(n log n) algorithm was given by Preparata and Shamos [59]. In our case, since the half-planes are presorted by their polar angles, it can be done in linear time. 3.4.2 Algorithm 2 In this section, we show how Chan’s output-sensitive convex hull algorithm ([15], Section 2.2) can be adapted to generate guaranteed hulls of aligned similar uncertain regions. Its running time is O(n log h), where n is the number of aligned similar uncertain regions, and h is the number of hull bitangents they induce. We are perhaps cheating a little, by considering the algorithm’s output to be the hull bitangents, and not the guaranteed hull that results from the intersection of their left half-planes. The difference between the complexities of these two types of output can be dramatic, as Figure 3.9 showed. In any case, the running time is an improvement over that of our first algorithm. Interestingly, Chan’s algorithm can also be used to generate the envelopes of functions in an output-sensitive manner ([15], Section 2.3). Since the hull bitangents of aligned similar uncertain regions can be represented as the vertices of an envelope of functions (by using the techniques of [55], Section 5.2, and [54], Section 4.2), this use of Chan’s algorithm will also construct guaranteed hulls in O(n log h) running time.4 In this section, we show how Chan’s algorithm can generate guaranteed hulls directly, in a more explicitly geometric fashion. Chan’s algorithm incorporates the Jarvis march algorithm with a novel grouping technique to construct the convex hull of a set of n points. In the first, or grouping step, the points are divided into dn/me groups of at most m points each. The convex hull of each group is calculated, using any O(m log m) convex hull algorithm (e.g., Graham’s scan [29]). The second, or wrapping step5 , is a variant of a Jarvis march. The hull vertices are generated in ccw order, by querying each group 4 Another output-sensitive envelope algorithm that can be used for this purpose is presented by Nielsen and Yvinec [56]. 5 The term ‘wrapping’ is chosen due to the step’s similarity to the gift wrapping algorithm, which in the planar case, is a Jarvis march. 93 for a candidate next vertex, then determining the actual next vertex from these candidates. These queries can be performed in constant time for each group (plus a term linear in the total size of the groups that is amortized over all the queries);6 hence the running time of the wrapping step is O(h · n m + n). If h were known in advance, then setting m to h would result in an O(n log h) algorithm. Since this is not the case, the algorithm is repeated with the following 1 2 sequence of ‘guesses’ at h: (22 , 22 , . . . , n). If the number of hull vertices generated for a particular guess does not exceed the guess itself, then clearly every hull vertex has been found, and the iterations can stop. It can be shown that the total running time of all of these iterations, and hence the total algorithm, is O(n log h). To describe how the algorithm can be modified for our purposes, we require some additional terminology. As before, we will assume that the regions are discs. Let Hθ represent the ordered sequence of upper hull bitangents of D whose polar angles do not exceed θ. Let {G1 , . . . , Gm } represent the partitioning of D into groups of dn/me discs. Each of these groups has a sequence of upper hull bitangents (B•e , . . . , Bw• ), ordered by polar angle. Note that each group has its own east and west discs, and each sequence starts and ends with their corresponding vertical hull tangents. It will be useful to think of a group’s sequence of bitangents as a stack, and define popping this stack as equivalent to removing its first (i.e., topmost) bitangent. The major work of the grouping step is constructing the guaranteed hulls of each group. Since we can use our first algorithm (Section 3.4.1) for this task, we will focus on adapting Chan’s wrapping step to our setting. We will start with Hθ = H π2 , and extend it by one bitangent at a time, until it is equal to H 3π . We 2 can construct H π2 by examining each group’s stack and retaining the rightmost bitangent B•e . Let us consider how to perform the remaining iterations efficiently. If, at the start of a wrapping step, the previous bitangent has the form B∗q , then the next bitangent in the sequence will (by Lemma 3.2) be of the form Bqv . We will show that not only will some group’s stack contain a bitangent of the form Bv∗ , but that if Bv∗ is not the topmost element of its stack, then we can recognize this fact and permanently remove those elements above Bv∗ . What is more, we can 6 The version of Chan’s algorithm that we discuss here incorporates one of his suggested improvements (‘Idea 5’; [15], p. 29). 94 do this without knowing the actual value of v. To elaborate, at the start of each wrapping step, we process each group’s stack. If it is nonempty, we pop bitangents from it while they satisfy at least one of the three conditions of the following lemma: Lemma 3.13 Let B∗q be the last bitangent added to Hθ , Gi be any group, and Bxy be the topmost bitangent on Gi ’s stack (with x 6= q). If (i) Dx contains Dq , (ii) θqx ≤ θ∗q , or (iii) θqx ≥ θxy , then either Dx does not contribute any upper hull bitangents of D with polar angle greater than θ∗q , or some bitangent Bxz exists lower in Gi ’s stack. Proof. Let B∗q , Gi , and Bxy be as described. If Dx contains Dq , then (by Lemma 3.11) Dx does not contribute any hull bitangents, and we are done. If θqx ≤ θ∗q , then θqx < θ∗q , and (by Lemma 3.4) θ∗q ∈ [θxq θqx ); hence θ∗q ≥ θxq > θqx . Now observe that by Lemma 3.4, no hull bitangent Bx∗ of D can exist with polar angle greater than θ∗q , and we are done. Suppose θqx ≥ θxy . If no hull bitangent Bx∗ of D exists where θx∗ > θ∗q , we are done. Otherwise, Bx∗ must be a hull tangent of Gi (since Dx ∈ Gi ). Lemma 3.2 then implies that there exists some hull bitangent Bxz of Gi , where θxz ≥ θx∗ . Lemma 3.9 now implies that θxz > θqx > θ∗q , and since θxy ≤ θqx , Bxz must be a distinct bitangent appearing below Bxy in Gi ’s stack, completing the proof. We now show that after each group’s stack has been processed in this way, the next hull bitangent of D can be deduced by examining the topmost elements of the remaining nonempty stacks. Suppose we are at start of an iteration of the wrapping step, that B∗q was the last hull bitangent added to Hθ , and that the next hull bitangent to be added is Bqv . For the moment, assume Dq and Dv belong to different groups. We will first demonstrate that some group’s stack contains a bitangent of the form Bvw , where θvw > θqv . Observe that since Dv supports a hull bitangent of D, it must support a hull tangent of the group Gi containing Dv . Hence, by Lemma 3.2, Gi ’s stack must have originally contained some Bvw where θvw > θqv . We need only show that Bvw is still in this stack. If this is not the case, then Bvw must have been popped 95 during some previous iteration of the wrapping step. Let B∗q0 be the hull bitangent added to Hθ just prior to Bvw being popped. By Lemma 3.13, q 0 6= v, and one of these three conditions must hold: (i) Dv contains Dq0 , (ii) θq0 v ≤ θ∗q0 , or (iii) θq0 v ≥ θvw . Now, (i) must be false, otherwise Dv could not support a hull tangent of D. Also, hull bitangents B∗q0 and Bqv are supported by discs Dq0 and Dv respectively, so by Lemma 3.9, θqv ≥ θq0 v > θ∗q0 ; and since θvw > θqv , both (ii) and (iii) must be false. Now that we have shown that a bitangent of the form Bvw remains in Gi ’s stack, we will show that it is this stack’s topmost element (or at least, the topmost element is of this form as well). Since Bqv is a hull bitangent of D, it must also be a pure hull tangent of Gi . Lemma 3.2 implies that Gi ’s stack originally contained some bitangent Bv∗ , where θv∗ > θqv ; and as we have shown, Bv∗ must still be in the stack. Assume, by way of contradiction, that the stack’s topmost element is Bxy , where x 6= v. Since Bxy was not popped in the previous iteration, θxy > θqx > θ∗q . Note that Dq can support no bitangents with polar angle in the range (θ∗q θqv ); hence θqx > θqv , which implies θxy > θqv . Let L be Bqv , and M the directed line parallel to L that has right-tangent Dx ; see Figure 3.12. Observe that (i) Dv is the only disc of Gi right-tangent to L, (ii) M is strictly to the left of L, and (iii) θ(M ) < θxy . If we now rotate L and M ccw in tandem around Dv and Dx respectively, then by the argument within the proof of Lemma 3.2, L must reach some hull bitangent Bvz of Gi before M reaches Bxy , since Dv lies strictly to the right of M (preventing it from being a hull bitangent of Gi ) until L and M coincide at Bvx . Hence, θxy > θvz > θqv > θ∗q . As we showed above, this implies that Bvz must still be on Gi ’s stack, and Bxy is not its topmost element (a contradiction). Now consider the case where Dq and Dv both belong to Gi . Since B∗q was the last hull bitangent generated, Gi ’s stack must have a topmost element of the form Bqp . Suppose, by way of contradiction, that p 6= v. Since Dq was chosen to support the last hull bitangent, Bqp must not have been popped per Lemma 3.13. It follows that θqp > θ∗q . Note also that since B∗q and Bqv are consecutive hull bitangents of D, Dq can support no bitangent with polar angle within (θ∗q θqv ). Hence, θqp > θqv . In addition, since {Dq , Dv } ⊆ Gi , Dp can support no hull bitangent of Gi with polar angle in (θ∗q 96 θv∗ ); thus θqp > θv∗ , and we have a L = Bqv Dq D∗ Bq∗ Dv Dx M Si Dy Bxy Figure 3.12: Bxy cannot be the stack’s topmost element. contradiction (since this implies that θqp is not the topmost element). To conclude, we have shown that after popping each group’s stack in accordance with Lemma 3.13, the topmost element of Dv ’s group’s stack will be either Bqv or Bv∗ . We can thus examine the topmost element of each (nonempty) stack and generate a candidate bitangent Bqx . The candidate with the smallest polar angle must be Bqv , the next hull bitangent in Hθ . Each iteration of the wrapping step finds the next hull bitangent in time linear in the number of groups, plus the time spent popping each group’s stack per Lemma 3.13. Since a bitangent can be popped only once, the total time spent in the wrapping step is linear in n. It follows that the running time of the complete algorithm is O(n log h). 3.4.3 Guaranteed Hulls in R3 While the focus of this thesis is on planar regions of imprecision, it is worth describing how some planar guaranteed hull algorithms can be adapted to construct guaranteed hulls in 3-dimensional space (R3 ). In particular, we will provide algorithms for constructing guaranteed hulls for uncertain regions that are (i) translations of a common base shape, (ii) spheres (i.e., the R3 -analogs of discs), and (iii) axis-aligned boxes (the R3 -analogs of rectangles). To work with guaranteed hulls in R3 , we need some terminology. A plane P 97 in R3 can be specified by a point and a normal vector. We will refer to the closed half-space bounded by P that does not contain the normal vector as the lower halfspace of P , and denote it by P − . If a lower half-space P − intersects every region in a set D, and three of these regions do not penetrate P − , then we will say that P is a hull plane of D, and we will say that each non-penetrating region is upper-tangent to P (we can also say that P is upper-tangent to these regions). thus R3 -analogs Hull planes are of hull bitangents. They also have a fundamental relationship to guaranteed hulls: Theorem 3.14 If GH(D) has nonzero volume, then GH(D) is the intersection of the lower half-spaces of the hull planes of D. Proof. We will assume GH has nonzero volume. Let P 0 be the set of planes whose lower half-spaces intersect every region of D. We start by proving that GH = I(P 0 ) (where I(P 0 ) denotes the intersection of the lower half-spaces of P 0 ). If some point p ∈ GH is not within I(P 0 ), then p lies strictly above some plane P ∈ P 0 . Since every region of D intersects P − , there exists a feasible set S of D that lies below P . This is a contradiction, since it implies that p ∈ / CH(S). Hence, p ∈ I(P 0 ), and GH ⊆ I(P 0 ). To show that I(P 0 ) ⊆ GH, let p be any point in I(P 0 ), and suppose p ∈ / GH. Then for some feasible set S, p ∈ / CH(S), which implies that S lies strictly below some plane Q containing p. Let P be the plane of P 0 that is parallel to Q and whose lower half-space contains no other plane of P 0 . Observe that some Di ∈ D is upper-tangent to P . Since p ∈ I(P 0 ), P must lie on or above Q; but then S lies strictly below P , a contradiction (since no point of Di lies below P ). We complete the proof by showing that I(P 0 ) = I(P). We will do this by demonstrating that every boundary point s of I(P 0 ) lies on a plane of P. It then follows that P is nonempty, and (since I(P 0 ) and I(P) are convex, and I(P 0 ) ⊆ I(P)) I(P) = I(P 0 ) = GH. Let s be any boundary point of GH, and q any point interior to GH. There must exist a plane of P 0 that contains s; let P be the closest such plane to q, and t(P ) be the number of regions of D upper-tangent to P . If t(P ) = 3, then P is a hull plane, and we are done. Otherwise, P contains a line L through each of the (at most two) tangency points of P . Now observe that since q ∈ / L, we can pivot P 98 around L to reach a new plane of P 0 that is closer to q, or whose lower half-space excludes s. Since both cases imply a contradiction, we are done. As we observed in Section 3.3, for uncertain regions that consist of translated versions of a common shape, we can generate their hull bitangents by constructing the convex hull of a set of similar points drawn from the regions. The same approach applies in R3 . The convex hull of a set of similar points is a polyhedron, and can be found in O(n log n) time. The hull planes of the original regions can be found by first constructing the planes containing the faces of this polyhedron, then lowering each of these planes by an appropriate distance. For regions such as spheres and polyhedra with a small number of faces, this poses no challenge. The intersection of the hull planes derived from the convex hull can be found in an additional O(n log n) time [59]. The running time of the complete algorithm is thus optimal. The second R3 guaranteed hull algorithm works with uncertain spheres of varying sizes. We will first describe the algorithm of Rappaport ([61], Section 3) for constructing the guaranteed hull of uncertain discs, then show how it can be adapted to R3 . The algorithm applies a simple transformation to each disc, then uses any standard algorithm to construct the convex hull of the transformed discs. Next, it transforms each convex hull edge to a hull bitangent of the original discs. Let rmax be the radius of the largest disc in D, a set of uncertain discs. A transformed set of discs is constructed, D 0 = {D10 , . . . , Dn0 }, by replacing each disc’s radius, ri , with the transformed radius ri0 = rmax − ri . It is easy to show that the hull bitangents of the original discs correspond to the boundary segments7 of the convex hull of these transformed discs. This algorithm extends directly to R3 . We transform the original spheres, construct their convex hull, lower each plane that contains a face of this hull by distance rmax , and output each lowered plane as a hull plane. By Theorem 3.14, the intersection of the half-spaces of these hull planes yields the guaranteed hull of the original spheres. The convex hull of a set of spheres can be constructed in O(n2 ) time [11]. 7 The other components of the hull’s boundary, arcs on the boundaries of the transformed discs, are ignored. 99 While the boundary of this hull is composed of planar facets, conical facets, and spherical facets [11], we are interested only in the planar facets. The running time of the complete algorithm is thus O(n2 +m log m), where m, which never exceeds n, is the number of planar facets. Our third R3 algorithm, which constructs the guaranteed hull of uncertain axisaligned boxes (with varying sizes in each dimension), is an adaptation of Nagai et al’s algorithm for guaranteed hulls of uncertain axis-aligned rectangles [55]. Their algorithm is based on the observation that any hull bitangent of a set of axis-aligned rectangles will be tangent at points that correspond to the same corners of two of the rectangles. The algorithm calculates the convex hulls of four sets of points, where each set consists of a particular corner drawn from each rectangle (e.g., the upper left corners), then constructs the guaranteed hull as the intersection of these four convex hulls. Converting this algorithm to R3 is simply a matter of constructing the intersection of eight convex hulls, corresponding to the eight corners of the axis-aligned boxes. The eight convex hulls can be constructed in O(n log n) time, and the intersection of the lower half-spaces corresponding to these convex hulls’ faces can be constructed in O(n log n) time as well, yielding an optimal algorithm. 3.5 Possible Hull The possible hull of a set8 of uncertain regions D is defined formally as the set of points p such that for some feasible set S of D, p ∈ CH(S). We will denote the possible hull of D by P H(D) (or, when D is clear from the context, simply by P H). One use for possible hulls is in identifying which regions of D are guaranteed extreme regions: regions that contribute a vertex to the convex hull of every feasible set of D. Theorem 3.15 If D is a set of at least three uncertain regions, then Di ∈ D is a guaranteed extreme region of D iff it does not intersect the possible hull of D \ {Di }. Proof. Let Di be a region of D, where |D| ≥ 3. If Di does not intersect P H(D \ 8 Throughout this section, we will assume that sets are multisets. 100 {Di }), then Di intersects the convex hull of no feasible set of D \ {Di }. Since this implies that Di is linearly separable from every such feasible set, Di must be a guaranteed extreme region of D. Suppose Di intersects P H(D \ {Di }). If there exists a feasible set S of D \ {Di } such that Di intersects CH(S) at a point p that is not a vertex of CH(S), then p is not a vertex of CH({p} ∪ S), Di is not a guaranteed extreme region of D, and we are done. If no such S exists, then Di must intersect a vertex p of the convex hull of some feasible set S 0 of D \ {Di }. Let Dj be the region of D \ {Di } contributing vertex p to CH(S 0 ). Now, Di must penetrate Dj , otherwise (since p ∈ Di ∩ Dj ) Di and Dj violate the regions’ valid intersection property. It follows that some point q ∈ Di lies in the interior of a line segment connecting Dj and some third region of D. Since q is not a vertex of the convex hull of any feasible set which includes this segment’s endpoints, Di is not a guaranteed extreme region of D. Nagai et al [55] have observed that for convex uncertain regions, the convex hull of the regions and their possible hull are the same. Hence, the guaranteed extreme regions of a set of convex uncertain regions are exactly those regions that lie outside the convex hull of the other regions in the set. In fact, it can be shown that in this case, a region is a guaranteed extreme region iff it lies outside the convex hull of its Voronoi neighbors. This implies that for many natural types of uncertain regions, such as discs and convex polygons with a small constant number of vertices, the guaranteed extreme regions can be identified in O(n log n) time. In this section, we develop a sequence of possible hull algorithms. Each algorithm in the sequence will use its predecessor as a subroutine, and the final algorithm will construct the possible hull of a set of (possibly nonconvex) polygons. We begin by stating some properties of possible hulls. Lemma 3.16 P H({A}) = A. Lemma 3.17 P H({A, B}) = [ ab. a∈A,b∈B Lemma 3.18 [ Di ⊆ P H(D). Di ∈D 101 Lemma 3.19 If A and B are nonempty sets of uncertain regions, then P H(A ∪ B) = P H({P H(A), P H(B)}). Proof. Let Q = P H({P H(A), P H(B)}). Suppose p is a point within P H(A ∪ B). Then there exists a feasible set S of A ∪ B such that p ∈ CH(S). Let Sa and Sb be the subsets of S corresponding to the subsets A and B. By using the definition of convex hull, it is easy to show that (i) CH(S) = CH(CH(Sa ) ∪ CH(Sb )), and (ii) there exist points a and b within CH(Sa ) ∪ CH(Sb ) such that p ∈ ab. Now, without loss of generality, either (i) a, b ∈ CH(Sa ), or (ii) a ∈ CH(Sa ) and b ∈ CH(Sb ). If (i), then ab ⊆ CH(Sa ), and since (by definition) CH(Sa ) ⊆ P H(A), ab ⊆ P H(A), which implies (by Lemma 3.18) that ab ⊆ Q. If (ii), then since CH(Sa ) ⊆ P H(A) and CH(Sb ) ⊆ P H(B), Lemma 3.17 implies ab ∈ Q. Hence P H(A ∪ B) ⊆ Q. If p is a point in Q, then by Lemma 3.17, p ∈ ab, where a ∈ P H(A) and b ∈ P H(B). There must then exist feasible sets Sa of A and Sb of B where a ∈ CH(Sa ) and b ∈ CH(Sb ). Note that ab is within CH(Sa ∪ Sb ), and since Sa ∪ Sb is a feasible set of A ∪ B, CH(Sa ∪ Sb ) is within P H(A ∪ B); hence Q ⊆ P H(A ∪ B). Lemma 3.20 The possible hull of any set of two or more connected uncertain regions is simply connected. Proof. Let D = {A, B} be a set of connected uncertain regions (Lemma 3.19 implies that the proof extends by induction for sets of more than two regions). By Lemma 3.17, every point of P H is connected within P H to both A and B; and by Lemma 3.18, both A and B lie within P H. Hence P H is connected, and we need only show that it has no holes. We will do this by showing that every exterior point of P H is the source of a ray exterior to P H. Let q be any point exterior to P H. We first make the observation that if no lines through q that are tangent to A exist, then (i) every line through q will intersect A, and (ii) at least one line through q will intersect A to both sides of q. Since the same argument applies to B, if no such lines exist for A or B, then some segment ab exists (where a ∈ A and b ∈ B) that contains q, implying that q ∈ P H, a contradiction. Hence, without loss of generality we can assume that there exist 102 directed lines L1 and L2 through q that are right-tangent to A. Let W1 (resp., W2 ) be the wedge lying on or to the right (resp., left) of both L1 and L2 . Note that W1 contains A. Note also that W2 cannot intersect B, otherwise some segment ab exists that contains q. Let R be any ray from q lying in W2 . Observe that no point → r ∈ R can lie on a segment ab, since for any choice of a ∈ A, the portion of ray − ar lying at or beyond r lies within W2 (Figure 3.13). Hence, by Lemma 3.17, R does not intersect P H. A a W1 L2 B q R r W2 L1 Figure 3.13: Lemma 3.20 Corollary 3.21 Let (p1 , . . . , pk , p1 ) be a cyclic sequence of points. If each consecutive pair (pi , pi+1 ) are endpoints of a segment known to lie within P H, and P is a simple polygon whose edges lie on these segments, then the interior of P lies within P H. A region R is star-shaped if there exists a point s ∈ R such that ∀r∈R sr ⊆ R. The kernel of a region R is the set of all such points s that have this property. 103 Hence, a region is star-shaped if it has a nonempty kernel. Theorem 3.22 The possible hull of any set of two or more connected uncertain regions is star-shaped. Proof. Let D = {A, B} be a set of two connected uncertain regions (Lemma 3.19 can be applied to prove the claim for sets of more than two regions). We will prove that P H is star-shaped by showing that it has a nonempty kernel. Let I = CH(A) ∩ CH(B). If I 6= ∅, then let s be any point in I. Observe that there must exist points a1 , a2 ∈ A, and b1 , b2 ∈ B where s lies on both a1 a2 and b1 b2 . Let q be any point in P H. By Lemma 3.17, there exist points a0 ∈ A, b0 ∈ B such that q ∈ a0 b0 . Without loss of generality we can assume that a1 , b1 , → and a0 are on or to the left of − sq, and that a , b , and b0 are on or to the right of 2 2 − → sq. There are now two cases (Figure 3.14): s, b1 , and a2 are either to the left or −−→ to the right of a1 b2 . In both cases, we can construct a cyclic sequence of points defining a polygon (per Corollary 3.21) that lies within P H, and which contains sq. (In the former case, the sequence is (a1 , b2 , a2 , b0 , a0 , b1 ); and in the latter case, it is (b1 , a2 , b0 , a0 ).) Since this holds for any q ∈ P H, s is in the kernel of P H. a1 a0 b1 a0 b1 a1 s b2 (i) q a2 s b0 q a2 (ii) b2 b0 Figure 3.14: Theorem 3.22 If I = ∅, then there exists line L1 right-tangent to A and left-tangent to B, and line L2 left-tangent to A and right-tangent to B. Let s be the point where L1 and L2 cross, and q be any point in P H. By Lemma 3.17, q ∈ ab, for some a ∈ A, b ∈ B. Note that there exist points a0 ∈ A and b0 ∈ B such that s ∈ a0 b and s ∈ ab0 . We 104 can again construct a sequence of points that defines a polygon that lies within P H and contains sq. If s lies to the right of ab, this sequence is (b, a, b0 , a0 ); otherwise, it is (a, b, a0 , b0 ). From this point on, we will assume that each of the uncertain regions is a polygon (or a point, which can be viewed as a degenerate polygon). The possible hull of one such set is shown in Figure 3.15. Note that as a consequence of our general position assumption, no three vertices lie on the same line, which implies that any pair of polygon edges will intersect only at their endpoints (if the edges belong to the same polygon) or at a single interior point (if the edges belong to distinct polygons). Figure 3.15: Possible hull of uncertain polygons Before proving a bound on the complexity of possible hulls of polygons, we will require some terminology. We will refer to an edge of a polygon of D as a native segment, and to a segment connecting vertices of distinct polygons of D as a bridge segment. Lemma 3.23 Each point on the boundary of P H lies on either a native segment or a bridge segment of D. Proof. Let q be a point on the boundary of P H. Note that q must lie on the boundary of some feasible hull C of D, and cannot lie in the interior of any feasible 105 hull of D. Suppose, by way of contradiction, that q does not lie on a native or bridge segment of D. If q is a vertex of C, then q is some feasible point within some polygon A ∈ D. Since q does not lie on a native segment of D, it must lie in the interior of A. Observe that we can replace q by some nearby q 0 ∈ A to create a new feasible hull with q in its interior. Since P H contains this feasible hull, q cannot lie on the boundary of P H, a contradiction. Point q must therefore lie in the interior of some segment ab, where a and b are vertices of C, and a and b lie within distinct polygons A and B of D. We can assume that a (resp., b) is the farthest point of A (resp., B) that intersects the line containing ab. Note that this implies that a lies on some edge a1 a2 of A, and that b lies on some edge b1 b2 of B. Since q is not on a bridge segment, we can assume (without loss of generality) that a is not a vertex of A, and that neither a1 nor a2 are collinear with ab. Now observe that the sequence (a1 , a2 , b) defines a polygon (per Corollary 3.21) within P H whose interior contains q, a contradiction. Theorem 3.22 and Lemma 3.23 imply that P H is a star-shaped polygon. As we now show, it has no more vertices than D. Theorem 3.24 If D is a set of uncertain polygons with a total of n vertices, then the boundary of P H(D) has at most n vertices. Proof. Let D be a set of uncertain polygons. We will show that each vertex v on the boundary of P H that is not already a vertex of D can be associated with a subset of the boundary of one of the polygons A of D that (i) is disjoint from the subset associated with any other vertex of P H, and (ii) contains a vertex of A (in the interior of P H) to which we can charge v. If v is a convex vertex, then since the native and bridge segments of D lie within P H (Corollary 3.21) and contain its boundary points (Lemma 3.23), v must be a vertex of D, and we are done. If v is a reflex vertex, we proceed as follows. Let u (resp., w) be the vertex of P H immediately preceding (resp., following) v on the boundary of P H. Let s1 s2 (resp., t1 t2 ) be the native or bridge segment containing uv (resp., vw). If s1 s2 and t1 t2 are native segments of the same polygon of D, then v(= s2 = t1 ) must be a vertex of this polygon, and we are done. If s1 s2 and t1 t2 106 are native segments of distinct polygons of D, then sequence (t1 , s1 , t2 , s2 ) induces (by Corollary 3.21) a polygon within P H whose interior contains v, implying that v is not a vertex of P H, a contradiction. Without loss of generality, there are two remaining cases: s1 s2 is a native segment of a polygon A ∈ D, and t1 t2 is a bridge segment; or s1 s2 and t1 t2 are both bridge segments. Consider the former case. We will let A represent the set D \ {A}. Note that no point p ∈ A can lie to the right of s1 s2 , otherwise (by Corollary 3.21) triangle 4ps2 s1 is within P H, and v is interior to P H. Hence, t2 ∈ A, and t1 ∈ A. Assume, for the sake of contradiction, that w 6= t2 . Then the next boundary edge of P H is wx, lying on some native or bridge segment q1 q2 , where q1 is to the left of t1 t2 and q2 is to its right (Figure 3.16). If q2 ∈ A, then 4wt1 q2 is within P H, and v is interior to P H. Thus q2 ∈ A, and lies to the left of s1 s2 . But now segments wv, vs1 , and q2 w bound a triangle within P H whose interior contains v, and we have a contradiction. Hence, w = t2 . Now note that the portion of the boundary of A from v ccw to w exclusive lies in P H’s interior, and includes vertex s2 of A. We can thus charge v to s2 . q2 x t2 w v s1 u PH s2 A t1 q1 Figure 3.16: Theorem 3.24 Consider the latter case: s1 s2 and t1 t2 are both bridge segments, and s1 is a vertex of a polygon A ∈ D. We can use very similar arguments to that of the first case to show that u = s1 , and w = t2 ∈ A. There is now a path on the boundary of A, ccw from u to w exclusive, that is interior to P H and contains some vertex of A to which we can charge v. 107 It is easy to construct examples of polygons for which the bound of Theorem 3.24 is tight. See, for example, Figure 3.17. Figure 3.17: Uncertain polygons with n vertices generating P H with n vertices 3.5.1 Point and Polygon The first of our possible hull algorithms constructs the hull of a point and a polygon. Suppose P is a (possibly nonconvex) simple polygon, and s is a point. We wish to construct P H = P H({s, P }); see Figure 3.18. We will maintain our general position assumption that s and the vertices of P are distinct, and that no three of them are collinear. P s Figure 3.18: Possible hull of polygon P and point s Our algorithm will actually construct the possible hull of a simple polygonal chain and a point, which is sufficient for two reasons: by Lemma 3.20, the possible hull of P and s is equal to the possible hull of the boundary of P and s; and the boundary of P is a simple polygonal chain (albeit one with no endpoints). 108 The algorithm may be useful in other contexts as well, since the possible hull of a polygon P and a point s is the smallest polygon containing P whose kernel contains s. Our algorithm is reminiscent of Melkman’s algorithm for finding the convex hull of a simple polygonal chain [50] for the following reasons. First, both are on-line algorithms: they require no preprocessing of the input, and maintain the respective hulls (i.e., convex vs. possible) of the portion of the input seen so far. Second, both algorithms look for points where the polygonal chain enters and emerges from the interior of the hull, and rely upon the simplicity of the polygonal chain to perform this efficiently. We motivate our algorithm with the following observation, implied by Lemma 3.17: the possible hull of a point s and a polygonal chain P is equal to a union of triangles, where each triangle’s vertices are s and the endpoints of an edge of P . Let (p1 , . . . , pn ) be the ordered vertices of P . We will denote the connected subset of P from a to b by (a . . . b). Our algorithm starts with point u initialized to p2 , and H initialized to the possible hull of s and (p1 . . . p2 ), which (as noted above) is the triangle with vertices s, p1 , and p2 . The algorithm advances u along P , processing each new edge (or part of an edge) in one of two ways. If the edge is exterior to H, then we expand H by adding its associated triangle; and if the edge is interior to H, then we skip the edge and advance u until it emerges from H. In either case, at the start of each iteration, u is a point that is on both P and the boundary of H. When u reaches pn , H will equal P H({s, P }). Since H is a simple polygon, we will assume its vertices have a ccw ordering. For added flexibility when manipulating H, we will associate with it a variable orientation whose values are ccw or cw. If a and b are adjacent vertices of H, with b ccw from a, then we will consider b to follow a when H has ccw orientation, and to precede a when H has cw orientation. Our algorithm consists of these steps: 1. Initialize H to be the triangle with vertices s, p1 , and p2 , and set u to p2 . 2. If u = pn , stop. 3. Set the orientation of H to ccw (resp., cw) if s lies to the left (resp., right) of 109 uv, where v is the vertex of P following u.9 4. If the points of P immediately following u lie in the interior of H, go to step 6. 5. Expansion step. Let T = 4suv, uv’s contribution to the hull. Starting with x = u, advance x along H, deleting those edges that lie within T , until the first of three events occurs: (i) H has no more edges (Figure 3.19). Replace H with T , advance u to v, and go to step 2. v v T u p1 pn T u p1 pn H H s s Figure 3.19: Expansion step, case (i), before and after changes to hull. (ii) x reaches the point where edge ab of H intersects sv. Replace ab with edges uv, vx, and xb (Figure 3.20). Advance u to v, and go to step 2. pn v pn T p1 v T u p1 b x a H H s u s Figure 3.20: Expansion step, case (ii), before and after changes to H. 9 In each of the figures that follow, H has ccw orientation according to this rule. 110 (iii) x reaches the point where edge ab of H intersects uv. Replace ab with ux and xb (Figure 3.21). Advance u to x, and go to step 2. b b T v v u pn pn a p1 T u x a H p1 H s s Figure 3.21: Expansion step, case (iii), before and after changes to H. 6. Interior step. Starting with x = u, advance x along P until the first of two events occurs: (i) x reaches pn ; stop. (ii) x reaches the point where P emerges from H’s interior (Figure 3.22). Split cd, the edge of H containing this point, into edges cx and xd. Advance u to x, and go to step 2. H p1 d v u x pn s, c Figure 3.22: Interior step, case (ii). We now examine the steps of the algorithm to prove that it generates the possible hull of s and P . We will use induction to show that at the start of each iteration of the algorithm, the following invariants hold: 111 1. H is the possible hull of s and (p1 . . . u). 2. u lies on the boundary of H. The invariants clearly hold for the base case, since the initial hull, H, is triangle 4sp1 u (where u = p2 ), and is thus equal to P H({s, (p1 . . . u)}). Suppose the invariants hold for every iteration until u reaches a particular position along P . We first consider the case where an expansion step is to be performed. If every edge of H lies within T (case i), then T is the possible hull of (p1 . . . v), and the invariants are satisfied. Suppose instead that some edge of H does not lie within T . Consider the the boundary points of H following u. Since s lies in the kernel of H, the first such point where the boundary crosses an edge of T must lie on sv (case ii), or in the interior of uv (case iii). In the former case, modifying the boundary of H as stated has the effect of expanding H to include 4suv(= T ); and in the latter, it has the effect of expanding H to include 4sux(⊂ T ). Since u is advanced to v in the former and x in the latter, we are thus adding exactly that portion of the hull contributed by those points of P between the old and new u, which satisfies invariant (1); and since the new u is not interior to H, invariant (2) is also satisfied. Let us now consider the case where an interior step is to be performed. We first show that we can ignore the points on (u . . . x). Since H is the possible hull of s and (p1 . . . u), and (u . . . x) ⊂ H, we have (by Lemma 3.17): P H({s, (p1 . . . x)}) = {sp | p ∈ (p1 . . . x)} = {sp | p ∈ (p1 . . . u)} ∪ {sp | p ∈ (u . . . x)} = H ∪ {sp | p ∈ (u . . . x)} =H. The only change we make to H is to split edge cd at x. As this does not actually change the boundary of H, and x (the new location of u) lies on this boundary, both invariants are satisfied. We have succeeded in proving the correctness of the algorithm. Before analyzing its running time, we will disclose some additional implementation details. Since P is immutable, we can store its vertices in an array. The vertices of H are 112 stored in a doubly-linked list, so that inserting or removing a vertex, or accessing a vertex’s neighbor, can be done in constant time. Since u occurs in both P and on the boundary of H, we maintain a pointer to its location in both data structures. In step 4, we need to determine if the points of P immediately following u lie in the interior of H. This is true iff the vertex of H following u is to the right (or, if the hull has cw orientation, to the left) of uv, and hence can be determined in constant time. In the interior step, we must determine which edge of H contains the point x where P emerges from H’s interior. To do this efficiently, we first characterize each boundary edge of H as being either a polygonal edge (one lying on an edge of P ) or a radial edge (one lying on a ray from s through a vertex of P ). Lemma 3.25 At the start of any iteration in which an interior step occurs, the following conditions hold: (i) exactly one of the edges of H incident with u is a radial edge; and (ii) if this radial edge is not incident with s, then x (the point where P emerges from H’s interior) must lie on this edge as well; otherwise, x must lie on an edge incident with s. Proof. At the start of an interior step, the edges of H incident with u cannot both be polygonal edges, otherwise (since an interior step is about to occur) this would imply that u is incident to three edges of P , which is impossible. Suppose instead that that they are both radial edges. Note that each radial edge of H is induced by a distinct vertex of P (unless a radial edge has just been split in step 6(d); but each such step is immediately followed by an expansion step that removes one of these two edges). Hence, s and two distinct vertices of P must be collinear, contradicting our general position assumption. To prove (ii), we first note that since P is simple, x must lie in the interior of cd, a radial edge of H. Let ab be the radial edge of H incident with u. Assume by way of contradiction that ab and cd are distinct edges, and that at least one of them is not incident with s. This edge must then be adjacent to (distinct) polygonal edges y1 and y2 (see, for example, edge ab in Figure 3.23). Now observe that (u . . . x) partitions H into two pieces, and since ab 6= cd, y1 and y2 must lie on opposite sides of (u . . . x). We now have a contradiction, since the interiors of 113 paths (p1 . . . u) and (u . . . x) must intersect, which implies that P is nonsimple. y1 H a c x b(= u) y2 d Figure 3.23: Lemma 3.25. Lemma 3.25 implies that in order to find x while moving along edges of P during an interior step, we need check for intersections of P with at most two edges of H: the single radial edge incident with the point of entry u, and (if that edge is also incident with s) the other edge incident with s. Every step of the algorithm can be done in constant time, except for the expansion and interior steps, which can be done in time proportional to the number of vertices that are: (i) visited on the chain; (ii) inserted into the hull; or (iii) removed from the hull. Since a vertex can be removed from the hull only once, the total running time of the algorithm is thus bounded by the number of iterations (which is O(n), since each iteration advances u to a distinct vertex of either H or P ), plus the number of vertices processed by the algorithm. The vertices processed include those of P , plus any vertices that ever appear in H. Consider the start of a particular iteration, where H is the current hull, u is the current position on P , and v is the vertex of P following u. We will prove a linear bound on the total vertices by showing that at most three vertices are introduced to H by the addition of triangle T = 4suv. Each new vertex (other than v) is a point where the boundaries of H and T cross, and hence must lie on either uv or vs (since su ⊂ H). Suppose for the sake of contradiction that two new vertices, p and q, lie in the interior of uv. We can assume that uv first enters H at p, then exits H at q. Since uv ⊂ P , and P is simple, → p and q must lie on radial edges of H. These radial edges must lie on rays − sa and 114 − → sb respectively, where a and b are vertices of P preceding u (Figure 3.24). The b v a H q p u P T s Figure 3.24: The interior of uv cannot contain two new vertices of H. → → or − path (a . . . b . . . u) (or (b . . . a . . . u)) in P cannot cross rays − pa qb (otherwise p or q would lie in H’s interior), nor can it cross pq (since P is simple). We now have a contradiction, since this implies that (a . . . b) is not connected to (u . . . v) within P . Hence, at most one new vertex lies in the interior of uv; and since s is in the kernel of H, at most one of the new vertices lies in the interior of vs. We can therefore conclude: Theorem 3.26 The above algorithm generates the possible hull of a point and a simple polygonal chain of n vertices in O(n) time. 3.5.2 Two Polygons Our second algorithm constructs the possible hull of a pair of uncertain polygons A and B. It starts with a polygon H equal to the convex hull of the pair, then modifies the boundary of H until H is equal to the possible hull of the pair. Such an approach seems promising, since CH(A, B) and P H({A, B}) will typically have many boundary edges in common. Consider Figure 3.25. If a boundary edge of CH(A, B) is a polygonal segment or a bridge segment, then (by Lemmas 3.17 and 3.18) it lies on the boundary of P H({A, B}) as well. Otherwise, its vertices must be nonadjacent vertices from the same polygon, and it will not be a boundary 115 aj ai A H B Figure 3.25: Initial polygon H; pocket lids shown as dashed lines. edge of P H({A, B}). As we manipulate H, both ai and aj will remain vertices of H, but the path (ai . . . aj ) on H’s boundary will change. We will refer to this path as a pocket of H, and to segment ai aj as the pocket’s lid. Our algorithm has these steps: 1. Initialize H to CH(A ∪ B). 2. Hull Contraction. For each pocket lid ai aj of H, perform the following steps: (a) Construct J, the possible hull of (ai . . . aj ) (on the boundary of A) and any point s from B. (b) Replace ai aj with path (ai . . . aj ) on the boundary of J. 3. Repeat the hull contraction steps with the roles of A and B reversed. 4. Hull Expansion. Perform the following steps: (a) Set u to h1 , any vertex of CH(A ∪ B). (b) Determine ray Tu as follows. If u is a vertex of A, then set Tu to the ray from u that is left-tangent to B; otherwise, set Tu to the ray from u that is left-tangent to A. 116 (c) If the vertex v of H following (i.e., in ccw direction) u is not left of Tu , then go to (f). (d) Let x be the point where the boundary of H next crosses Tu . (e) Replace boundary edges (u . . . x) of H with ux. (f) Advance u to the next convex vertex of H. If u 6= h1 , go to (b). 5. Repeat the hull expansion steps, substituting cw for ccw, and right for left. The hull contraction steps use the algorithm of the previous section to replace each pocket lid with a portion of a possible hull associated with the pocket. This contracts the hull by ‘taking bites’ out of the convex hull of the two polygons (Figure 3.26). Corollary 3.21 implies that after being modified in this way, each pocket lies within P H. As we will see, each pocket also has a particular monotonicity quality that both ensures the correctness of the remaining steps of the algorithm, and contributes to the algorithm’s efficiency. The hull expansion steps traverse the boundary of H, find tangent rays that potentially contain bridge segments of P H, and modify H to incorporate these segments (Figure 3.27). This has the effect of expanding the hull, by adding back some portions that were removed in the hull contraction steps. aj ai A s J H B Figure 3.26: Hull contraction: (ai . . . aj ) of J. pocket lid ai aj has been replaced by 117 u x tu Tu A h1 H B Figure 3.27: Hull expansion (ccw direction). We will prove that the algorithm is correct by showing that after all its steps have been performed, H is within P H, and P H is within H. To support these claims, we will need some additional terminology. We will denote by Tx the ray from point x that is left-tangent to a particular region (usually A or B), and by tx its point of tangency. If P = (u . . . v) is a path, and R is a region, then we will say that P is left-monotonic with R if (i) neither u nor v penetrate R, (ii) v lies in the right half-plane of Tu , and (iii) the sequence of rays from points of P that are left-tangent to R have nondecreasing polar angles.10 We define right-monotonicity analogously. One property of these paths is that any subpath can be replaced by a line segment, and monotonicity is preserved: Lemma 3.27 If path (u . . . v) is left-monotonic (resp., right-monotonic) with region R, then segment uv is left-monotonic (resp., right-monotonic) with R. Proof. We will prove only the case for left-monotonicity, since the two cases are symmetric. Suppose (u . . . v) is left-monotonic with R, and w is any interior point of uv. We need only show that the third property holds, since the truth of the first two does not depend upon w. Now, tv must lie in the left half-plane of Tw , 10 Points of P that lie on the boundary of R may lie on infinitely many lines tangent to CH(R). In this case, we can choose a tangent line that maintains the monotone property. 118 −→ −→ which implies that θ(wtv ) ≥ θ(Tw ); and since w ∈ uv, θ(Tv ) ≥ θ(wtv ). Thus, θ(Tv ) ≥ θ(Tw ). A symmetric argument shows that θ(Tw ) ≥ θ(Tu ). At the start of the hull expansion steps, the boundary of H consists of native segments, bridge segments, and pockets derived from the possible hulls J of the hull contraction steps. Together, these three sets of edges form a closed path that satisfies the conditions of Corollary 3.21; hence, at the start of step 4, H lies within P H. Observe that step 4 (and, by a symmetric argument, step 5) affects only pockets of H, for if uv is a native or bridge segment on the boundary of CH(A ∪ B), then v will not lie to the left of the ray constructed in step 4(b), and no change to H’s boundary will be made. Hence, in the rest of our analysis, we will focus mainly on the pockets of H. Lemma 3.28 At the end of the hull expansion steps, each pocket (ai . . . aj ) of H is both left- and right-monotone with B. Proof. We will show only that at the end of step 4, each pocket is left-monotone with B. A symmetric proof can then be used to show that step 5 induces rightmonotonicity, and (by Lemma 3.27) preserves left-monotonicity. The proof is by induction on the length of path (ai . . . u), a subpath of (ai . . . aj ). When u = ai , this path is a single point, exterior to B, and is thus left-monotone with B. Suppose (ai . . . u) is left-monotone with B, for some u ∈ (ai . . . aj ). The vertex v following u lies either to the left or right of Tu (Figure 3.28). Suppose v lies to the left of Tu . Observe that at the start of step 4, pocket (ai . . . aj ) is left-monotone with some point s ∈ B. Since the only changes made to the pocket by the hull expansion steps involve replacing subpaths with individual line segments, the pocket must (by Lemma 3.27) remain left-monotone with s throughout these steps. As a consequence, the pocket lies within triangle 4sai aj , and v lies to the right of us. Now, aj cannot lie to the left of Tu ; otherwise, θ(Tu ) < θ(ai aj ) < θ(Tai ), violating (ai . . . u)’s monotonicity with B. It follows that path (u . . . aj ) must cross Tu to reach aj , which implies the existence of point x of step 4(d). If x precedes tu on Tu , then θ(Tx ) = θ(Tu ); and if x follows tu on Tu , then we can rotate Tu ccw around x until it reaches Tx that is left-tangent to B, 119 aj aj ai ai v u Tu v tu B s u Tu Tv B s Figure 3.28: Lemma 3.28. which implies that θ(Tx ) ≥ θ(Tu ). In both cases, (u . . . x) is left-monotone with B. Now consider the case where v lies to the right of Tu . If we start with a line through v that is parallel to Tu , and rotate it ccw around v, it must become lefttangent to B before it has rotated by π radians; hence, θ(Tv ) > θ(Tu ), which implies that (u . . . v) is left-monotone with B. We can use a very similar argument to this second case to show that this monotonicity exists for every point following u (or x, if the first case applied) up to the pocket’s next convex vertex. Thus, after step 4(f), (u . . . v) is left-monotone with B. Let us prove our first claim, that H is within P H at the end of the hull expan- sion steps. We showed earlier that this is true at the start of these steps, so we will consider the changes made to a particular pocket (ai . . . aj ) in step 4 (the proof for step 5 is symmetric). Consider one iteration of step (e), the only step that changes H. Observe that u cannot have been added to H by a previous iteration of step (e), since u is convex and these steps can only add reflex vertices. Thus u must have been a vertex of J. Now, sai and saj are both within J, which implies that if s is a vertex of J it must be a reflex vertex. Hence, u 6= s. In fact, u must be a vertex of A. For if it is not, then by Lemma 3.23, u must be a point of intersection of L1 and L2 , where L1 and L2 are either edges of A, or bridge segments between vertices of A and s. Since u is neither s nor a vertex 120 of A, u lies in the interior of both L1 and L2 . But L1 and L2 both lie within P H, by Lemmas 3.17 and 3.18, which implies that u is not a convex vertex of J, a contradiction. In the proof of Lemma 3.28, we showed that point x of step (d) must exist. There are two cases (Figure 3.29): x is either strictly between u and tu , or it is at or beyond tu . We can now apply Corollary 3.21 to show that the region added to u x tu Ru u H H tu x B B Ru Figure 3.29: Expanding H remains within P H. H in step (e) lies within P H. In the former case, the region is bounded by path (u . . . x) and utu . In the latter case, since u and aj are both vertices of A, path (ai . . . aj ) on the boundary of A must cross Tu at a point y which is at or beyond x; hence, the region is bounded by segments utu , tu y, and path (u . . . x). We are now ready to prove our second claim, that P H ⊆ H. Let (ai . . . aj ) be any pocket of H, and let Li and Lj be the lines through ai and aj respectively that are left-tangent to B. These lines intersect at a point c to form a double wedge. B lies in one of these wedges, and Lemma 3.28 implies that pocket (ai . . . aj ) lies in the other wedge, within triangle 4cai aj (Figure 3.30). If some point x ∈ P H is not within H, then by Lemma 3.17, x lies on some segment ab, where a ∈ A, b ∈ B. Since x ∈ CH(A ∪ B), and all non-pocket edges of H lie on the boundary of CH(A ∪ B), x must lie within the triangle associated with some pocket; without loss of generality, we can assume that the pocket is (ai . . . aj ). Since B lies in the wedge opposite to 4cai aj , a must lie in this triangle. In fact, → − either a, or some point farther from b on ba, must be a point on pocket (ai . . . aj ). We can now apply Lemma 3.28 to show that x ∈ H, a contradiction. The algorithm is correct; we now show that it is efficient as well. Let n be the 121 aj Li A 4caiaj B c ai Lj Figure 3.30: Wedges induced by pocket (ai . . . aj ). total number of vertices of A and B. To perform step 1, we start by constructing the convex hulls of both A and B. This can be done in O(n) time, e.g., by using Melkman’s algorithm [50]. Once these hulls are known, the convex hull of A ∪ B can be constructed, in O(n) time, using the rotating calipers paradigm [72]. In the hull contraction steps, for each pocket ai aj , we construct the corresponding subset of the boundary from A, then calculate the possible hull of this boundary and an arbitrary point of B. Since each edge of A appears in only one of these subsets, and the possible hulls can be constructed in time linear in the size of the subset (Theorem 3.26), we can perform these steps in O(n) time. In step 4(b), we find the ray through a point u ∈ A that is left-tangent to B. Using reasoning very similar to that of the proof of Lemma 3.28, we can show that by remembering the vertices of CH(A) and CH(B) last used to calculate a ray tangent to one of these polygons, each iteration of this step can be performed in O(1 + m) time (by advancing the appropriate pointer until it reaches the vertex supporting a tangent line), where the total m over all iterations is linear in the complexity of the hulls, or n in the worst case. Since an edge can only be removed once, steps 4(d) and 4(e) can be performed in O(1 + k) time, where the total k over all iterations is linear in the number of vertices at the start of the expansion steps. By Theorem 3.24, H at this time (as well as at the start of step 5) has no more than n vertices; hence, the total k over all iterations is linear in n. Each of the remaining substeps of step 4 require constant time per iteration. Since each iteration advances u to a distinct vertex of H, there are O(n) iterations; thus the total time spent in 122 step 4 (as well as in step 5) is O(n). We have shown that each step of the algorithm can be performed in linear time. Hence, Theorem 3.29 The possible hull of a pair of polygons with n total vertices can be constructed in O(n) time. A natural question to ask is whether the possible hulls J need to be constructed within the hull contraction steps. It would be simpler, after all, if we could just contract each pocket ai aj to the corresponding path (ai . . . aj ) on the boundary of A. One reason why this cannot be done is that if this simplification is made, the pocket is not guaranteed to be left-monotone to any particular point of B at the start of the hull expansion steps, which impacts the efficiency of those steps. Consider the pair of polygons in Figure 3.31. Each of the ‘spikes’ in the boundary of A have a pair of convex vertices u that are not left-monotone with B. As noted earlier, the existing algorithm calculates Ru by advancing a pointer m times until it finds a vertex supporting tangent line Ru . If we use (ai . . . aj ) as the initial pocket, then m can be Θ(n) for each tangent line calculation, resulting in an O(n2 ) running time. We could improve the running time to O(n log n) by using a binary search to find Ru , but this is still inferior to the O(n) running time of the original algorithm. ai B A aj Figure 3.31: Lemma 3.28 will not hold. 123 3.5.3 Multiple Polygons We can leverage the results of the previous sections to produce an algorithm for the possible hull of any number of polygons, with minimal extra effort. Theorem 3.30 The possible hull of k polygons with a total of n vertices can be constructed in O(n log k) time. Proof. By Theorem 3.29, the possible hull for a pair of polygons can be constructed in O(n) time. Lemma 3.19 implies that we can apply this algorithm recursively, in a divide-and-conquer manner, to construct the possible hull of k polygons. In doing so, we increase the running time by a factor that is logarithmic in the height of a binary tree of k elements. 3.6 Closing Remarks This chapter was concerned with two variants of the convex hull problem involving uncertain sites: guaranteed hulls, and possible hulls. We provided an improved bound on the complexity of guaranteed hulls, and presented a novel algorithm for constructing guaranteed hulls of a particular class of uncertain regions. We also demonstrated how an existing output sensitive convex hull algorithm of point sites can be adapted for such regions, and showed how three planar guaranteed hull algorithms can be adapted to construct guaranteed hulls in R3 . In addition to identifying some properties of possible hulls of (not necessarily convex) uncertain regions, we presented an optimal algorithm for constructing possible hulls of uncertain polygons. To the best of our knowledge, these are the first results related to possible hulls of nonconvex uncertain regions. We can identify some open problems and areas of future research with respect to guaranteed hulls. All of the guaranteed hull algorithms that we have encountered (including the ones we have introduced here) operate by constructing a set of hull bitangents, then generating the guaranteed hull as the intersection of the bitangents’ half-planes. The complexity of this intersection may be much smaller than the number of bitangents (the intersection may even by empty), which suggests that any truly output-sensitive algorithm for guaranteed hulls will probably require a different approach. 124 Another open problem relates to the number of hull bitangents induced by a set of uncertain polygons. We have shown that the complexity of a guaranteed hull is strictly linear in the number of uncertain regions, regardless of the regions’ shapes (Theorem 3.6); however, the number of hull bitangents that produce this guaranteed hull is not necessarily linear. While we have proven linear bounds for the number of bitangents induced by aligned similar uncertain regions (Theorem 3.7), we have been unable to improve the O(nα(n)) bound given by [55] for uncertain polygons. The algorithms discussed in Section 3.4.3 rely on Theorem 3.14, and hence can only be relied upon to construct guaranteed hulls that are known to have nonzero volume. This situation could be improved by strengthening Theorem 3.14 to be a direct analog of Theorem 3.1 (we leave this as an open problem). In addition, the algorithm we presented for constructing guaranteed hulls of uncertain spheres is perhaps not optimal. Recall that one step of the algorithm involves constructing the convex hull of a set of spheres. As shown in [11], these hulls can have Θ(n2 ) complexity; but we are only interested in the planar facets of these hulls. We do not know if the number of these facets is always subquadratic. If so, then modifying the algorithm of [11] (or some other algorithm) to generate only planar facets of these hulls may yield a subquadratic guaranteed hull algorithm. We also leave as an open problem the complexity of guaranteed hulls of general uncertain regions in R3 . Figure 3.32: Applet for Chapter 3. An applet that implements some of the algorithms described in this chapter can be found at http://www.cs.ubc.ca/∼jpsember/uh.html. The applet (Figure 3.32) 125 includes the following sample files: 1. dischull.dat: Constructs the guaranteed hull of uncertain discs, using the algorithm of Section 3.4.1. 2. oshull.dat: Constructs the guaranteed hull of uncertain discs, using the algorithm of Section 3.4.2. 3. poly pt.dat: Constructs the possible hull of a point and a nonconvex polygon, using the algorithm of Section 3.5.1. 4. valley.dat: Constructs the possible hull of a pair of nonconvex polygons, using the algorithm of Section 3.5.2. 5. valley multiple.dat: Constructs the possible hull of several nonconvex polygons, using the algorithm of Section 3.5.3. 126 Chapter 4 Smallest Bounding Disc of Imprecise Points 4.1 Introduction In the 1-center problem (or, as it is sometimes called, the minmax location problem), we are given a set S of points in the plane, and we need to find the location of an additional point f such that the maximum distance from f to any point of S is minimized. The 1-center problem is a classic facility location problem from the field of operations research. In that setting, the points S are considered demand points, and the point f is considered the location of a facility to service the demand points. The 1-center problem has application in many other fields as well, such as computer graphics, mechanical engineering, biology, and even political science [28]. The 1-center problem can be generalized in many ways. For instance, the distance between the facility and demand points may use other metrics (e.g., the L1 or Manhattan metric), or may be replaced by a more complicated transportation cost, perhaps involving individual weighting factors on the demand points. Another generalization is the k-center problem. This differs from the 1-center problem in that more than one facility can be specified, and each demand point is free to use the nearest facility. In this chapter, we will consider only the 1-center problem, where the distances 127 are Euclidean (i.e., employ the L2 metric). In this setting, the 1-center problem is equivalent to what we will call the smallest bounding disc problem, in which the goal is to find the smallest disc that contains a point set. We will denote the smallest bounding disc of a point set S by B(S). There are two fundamental conditions that a disc C must satisfy in order to be a smallest bounding disc of S. The intersection condition holds if S is a subset of the disc, and the support condition holds if every section of the boundary of the disc that subtends an angle of at least π intersects the set. Clearly, the intersection condition must hold for C to be a bounding disc of S. The support condition must hold as well, for if some boundary arc of C subtends an angle of at least π and does not contain any points of S, then we can move the origin of C away from this arc and reduce its radius, yielding a smaller bounding disc of S. It is easy to show that any disc C that satisfies these conditions must be unique: if p is the origin of C, and q is any other point in the plane, then the support condition implies that at least one point of S must be strictly farther from q that it is from p. Hence, any disc centered at q that satisfies the intersection condition must be larger than C. An important property of the smallest bounding disc of a point set is that it is determined by a small number of points. If we start with a point set S, and repeatedly remove points from this set until no points can be removed without changing its smallest bounding disc, then we will be left with some subset of S. We will call this subset a critical subset of S. When it is clear from the context, we may refer to a member of a critical subset as a critical point. There may not be a unique critical subset for a point set. For instance, the vertices of a regular pentagon have five distinct critical subsets. We can make some observations concerning critical subsets. First, every critical subset of a point set S lies on the boundary of B(S). Second, every critical subset of S (where |S| ≥ 2) consists of either two or three points: if it has two, then the midpoint of the two critical points will be the 1-center of S (and the critical points are the endpoints of a diameter of the smallest bounding disc of S); and if it has three, then the angles subtended by the arcs between neighboring critical points will be less than π. The left set of points in Figure 4.1 has two critical points, and the right set of points has three. 128 Figure 4.1: Smallest bounding discs of point sets. This chapter is concerned with a variant of the 1-center problem, one involving imprecise points. We define the possible 1-centers of uncertain regions D to consist of those points in the plane that are the 1-center for some feasible set S of D. We will denote the possible 1-centers by C◦ (D) (or simply C◦ ). 4.1.1 Related Work The smallest bounding disc problem was first posed by Sylvester in 1857 [68]. A linear time algorithm for its solution was given by Megiddo [49], and Welzl provided a simple randomized algorithm to solve the problem in expected linear time [74]. Prior to Megiddo’s result, the fastest algorithm to solve this problem was given by Shamos and Hoey [65], who observed that the smallest bounding disc of a set of points must be centered on a vertex of the farthest Voronoi diagram of the points, which can be constructed in O(n log n) time. Some work has been done related to smallest bounding discs of uncertain points. The intersection radius of a set of uncertain regions is the radius of the smallest disc that is a feasible smallest bounding disc of the regions. Bhattacharya et al [10] showed that the intersection radius for a set of uncertain lines and points can be found in linear time. This result was later extended to sets of uncertain regions that may also include convex polygons, rays, and wedges [32]. For uncertain discs, Löffler and van Kreveld [48] have shown that an adaptation of Welzl’s algorithm [74] can find the intersection radius in expected linear time. The problem of computing the largest feasible smallest bounding disc for uncertain regions has been investigated by Löffler and Kreveld [48] as well. They 129 show that for uncertain discs or uncertain axis-aligned squares, the problem can be solved in linear time. 4.1.2 Contributions In this chapter, we introduce a new geometric figure: C◦ (D), the possible 1-centers of uncertain regions D. We identify some properties of this figure, and prove that its boundary points are contained by three families of curves induced by the uncertain regions. We demonstrate a connection between C◦ (D) and the guaranteed hull of D, one that can be exploited to optimize operations involving C◦ (D). For the cases where the uncertain regions D are discs or axis-aligned rectangles, we present two algorithms: an O(n log k) query algorithm to determine if a query point lies within C◦ , and an expected O((n + k 6 ) log k) construction algorithm to construct C◦ . The value k is the number of regions of D that are not interior to the guaranteed hull of D. This value never exceeds n, and can be substantially smaller than n. 4.2 Possible 1-Centers In this section, we identify some properties of the possible 1-centers of a set of general uncertain regions, and show that each boundary point of this figure lies on one of three types of curves induced by the regions. We also highlight a connection between possible 1-centers and guaranteed hulls. Finally, we look at the specific case where the uncertain regions are discs, and determine what form the three curves take in this case. Lemma 4.1 If D is a set of connected uncertain regions, then C◦ (D) is a connected, compact set. Proof. We can view the 1-center of a set of n uncertain regions as a mapping from a connected subset of R2n to R2 . This mapping is continuous, hence C◦ is compact ([70], Theorem 13.3.1) and connected (Ibid., Theorem 13.4.6). To describe C◦ , we will investigate where points on its boundary can lie. We will say that an element si of a feasible set S is a boundary site if it lies on the 130 boundary of its corresponding region Di . Lemma 4.2 Let D be a set of convex uncertain regions. If S is a feasible set of D whose 1-center lies on the boundary of C◦ (D), then two elements of S are boundary sites that also lie on the boundary of B(S). Proof. Let S be a feasible set whose 1-center, p, lies on the boundary of C◦ . We will first modify S by replacing any boundary sites si ∈ S lying in the interior of B(S) with points s0i interior to both B(S) and Di . Note that after this procedure, (i) every boundary site of S will lie on the boundary of B(S), (ii) B(S) (and p) have not changed, and (iii) all of the original boundary sites of S lying on the boundary of B(S) will remain. Hence, if we can prove the claim for the modified set, we have proven it for the original set. If S has no boundary sites, then observe that for any θ, there exists a nonzerolength vector with polar angle θ (denoted τθ ) such that each site si ∈ S can be replaced by the point si + τθ . This yields a new feasible set S 0 whose 1-center is p + τθ . As this is true for any choice of θ, p cannot lie on the boundary of C◦ ; a contradiction. If S has a single boundary site si , then si can be translated by some nonzerolength vector τθ , and each sj6=i ∈ S can be translated (if necessary) to an interior point of Dj so that the resulting set S 0 is (i) feasible, (ii) has the same 1-center as S, and (iii) has no boundary sites (Figure 4.2). We can now apply the previous argument to reach a contradiction. It is necessary at this point to define a variant form of tangency. We will say that region R is locally inside-tangent to disc C at point p if (i) p lies on the boundaries of both R and C, and (ii) every point of R within some distance > 0 of p lies within C. We may also say that C has local inside-tangent R at p. In Figure 4.3, polygon R is locally inside-tangent to disc C at p. Note that if R is a disc, then local inside tangency is the same as inside tangency. If Df and Dj are distinct regions of D, then we define the type I curve of Df and Dj to be the set of points that are the centers of discs which have outsidetangent Df and locally inside-tangent Dj . We denote this curve by If j . The subscript ‘f ’ is chosen since Df will usually be the farthest region from points of interest on this curve, as we will see later. 131 p si τθ s0i Di Figure 4.2: Lemma 4.2. C R p Figure 4.3: R is locally inside-tangent to C at p. Lemma 4.3 If D is a set of convex uncertain regions, and if p is both a point on the boundary of C◦ (D) and the 1-center of a feasible set S that has a critical subset of three points, then p must lie on a type I curve of D. Proof. Let p and S be such a point and feasible set respectively. We will show that p lies on a type I curve by proving that some region of D is outside-tangent to B(S), and another region of D is locally inside-tangent to B(S). If no region is outside-tangent to B(S), then observe that we can construct another feasible set S 0 with the same 1-center as S, by replacing each point of S on the boundary of B(S) with a point that is (i) slightly closer to p, and (ii) interior to its respective region. Lemma 4.2 then implies that p is not on the boundary of 132 C◦ , a contradiction. If no region is locally inside-tangent to B(S), then we can use a similar procedure to construct a feasible set S 0 of points that are slightly farther away from p, and again apply Lemma 4.2 to arrive at a contradiction. We now consider points s on the boundary of C◦ associated with feasible sets that have only two critical points. We will show that such points lie on one of two additional types of curves. We define the type II curve of distinct regions Di and Dj of D to be the boundary of the Minkowski average of the regions, and denote this curve by IIij . The Minkowski average of regions Di and Dj is the set 1 (pi + pj ) | pi ∈ Di , pj ∈ Dj . 2 We define the type III curve of distinct regions Df , Di , and Dj of D to consist of origins of discs Cp where (i) Cp has outside-tangent Df ; and (ii) p is the midpoint of boundary points of Di and Dj (see Figure 4.4). We denote this curve by IIIf ij . Dj Cp p Df Di Figure 4.4: Point p is on the type III curve of Df , Di , and Dj . Let D be a set of convex uncertain regions. We now show that each boundary point of C◦ (D) must lie on a type I, II, or III curve of D. Let p be a boundary point of C◦ that does not lie on a type I curve. There must exist a feasible set S whose 1-center is p; and by Lemma 4.3, every critical subset of S must have 133 exactly two points, with p as their midpoint. Let one such critical subset be Q = {qi ∈ Di , qj ∈ Dj }. We first consider the case where no third region Dk∈{i,j} of / D is outside-tangent B(S). Let us perform a procedure similar to that of the proof of Lemma 4.2, if necessary, to ensure that no sites other than qi and qj lie on the boundary of B(S). We can then apply Lemma 4.2 to prove that both qi and qj lie on their respective regions’ boundaries. Observe that there is a continuous, closed range of polar angles of directed lines through qi that are right-tangent Di . We denote this range by Right(i). Note also that each angle θ ∈ Right(i) has a corresponding angle θ + π, the polar angle of a directed line through qi that is left-tangent Di . We denote this corresponding range by Left(i). The Right(i) and Left(i) ranges are separated by disjoint, open ranges, representing the polar angles of lines through qi that are about to either leave Di , or enter the interior of Di , at that point. We denote the former range by Out(i), and the latter range by In(i). Note that for convex uncertain regions with smooth (that is, continuously differentiable) boundaries, both Right(i) and Left(i) contain a single value, while for other regions (e.g., polygons), this may not be the case; see Figure 4.5. Right(i) In(i) Right(i) In(i) Di Di I qi qi Out(i) Out(i) Left(i) Left(i) Figure 4.5: Polar angle ranges associated with qi ∈ Di . If Right(i) and Right(j) are not disjoint, then there exist parallel lines through qi and qj that are right-tangent Di and Dj respectively (Figure 4.6). It follows that p lies on the boundary of the Minkowski average of Di and Dj ; hence, p ∈ IIij . If Right(i) and Right(j) are disjoint, then it is easy to show that In(i) and Out(j) intersect. It follows that there exists an angle θ where θ ∈ In(i) and θ + π ∈ In(j) (Figure 4.7). Note that this implies that si and sj can be replaced by points 134 In(j) In(i) Right(i) Dj Di Right(j) qi Out(i) qj Left(j) Left(i) Out(j) Right(i) ∩ Right(j) Figure 4.6: Right(i) ∩ Right(j) 6= ∅. s0i , s0j interior to their respective regions, such that their midpoint remains at p. The resulting feasible set S 0 , whose 1-center is also at p, has no boundary sites. Lemma 4.2 then implies that p is not on the boundary of C◦ , a contradiction. Now suppose some Dk is outside-tangent B(S). To prove that p lies on IIIkij , We need only show that qi and qj are boundary sites. If this is not the case, then (without loss of generality) qi lies in the interior of Di . Note that Dj cannot be locally inside-tangent B(S), otherwise p lies on Ikj , which we know to be false. It follows that we can move qi and qj directly away from each other and into their respective regions’ interiors, to produce a feasible set S 0 , centered at p, but with no boundary sites. We can then apply Lemma 4.2 again to reach a contradiction. We have proven the following result. Theorem 4.4 If D is a set of convex uncertain regions, then every point on the boundary of C◦ (D) lies on a type I, II, or III curve of D. The guaranteed hull of a set of uncertain regions can be used to identify some regions that will not contribute to the possible 1-centers of the regions. 135 Right(i) In(j) In(i) Di Right(j) Dj qi qj Left(j) Out(i) Out(j) Left(i) In(i) ∩ Out(j) Figure 4.7: In(i) ∩ Out(j) 6= ∅. Theorem 4.5 If Di is a member of D, a set of uncertain regions, and Di lies in the interior of GH, the guaranteed hull of D, then C◦ (D) = C◦ (D \ {Di }). Proof. Let q be an interior point of GH. We start by showing that every directed line through q has at least one region of D lying strictly to its right. Let L be one such line. If no region of D lies strictly to the right of L, then there exists a feasible set S of D that lies in the left half-plane of L. But since GH ⊆ CH(S), q is not an interior point of GH, a contradiction. ◦ Let us denote C◦ (D \ {Di }) by C . Suppose, by way of contradiction, that ◦ some point p ∈ C is not within C◦ . Then there exists a disc Cp , centered at p, that is a feasible smallest bounding disc of D \ {Di } but not D. Note that this implies that Di does not intersect Cp , which means that there exists a directed line L that (i) is left-tangent to Di at a point q, and (ii) has Cp strictly to its left. Now note that since q is interior to GH, some region of D \ {Di } lies strictly to the right of L, and hence does not intersect Cp , a contradiction (since Cp must intersect every ◦ region of D \ {Di }). Hence, C ⊆ C◦ . ◦ Now suppose some point p ∈ C◦ is not within C . Then there exists a disc Cp , the smallest bounding disc of a feasible set S of D, that is not also the smallest 136 bounding disc of S \ {si }. This implies that Cp does not satisfy the support condition for S \ {si }, which implies that si is a boundary point of Cp . Now consider the line L left-tangent to Cp at si . Since si lies in the interior of GH, some region Dj ∈ D must lie strictly to the right of L. This is a contradiction, since it implies ◦ that sj ∈ S lies outside of Cp . Hence, C◦ ⊆ C . Figure 4.8: Possible 1-centers of uncertain discs. We now consider the possible 1-centers of uncertain discs (Figure 4.8). We start by showing how Theorem 4.5 can be applied to ‘filter out’ some of the discs that do not contribute any possible 1-centers: Lemma 4.6 Let D be a set of n uncertain discs. The k discs that do not lie in the interior of GH(D) can be found in O(n log k) time. Proof. We first use the algorithm of Section 3.4.2 to construct GH. The proof of Theorem 3.6 implies that h, the number of hull bitangents of D, is not more than k; hence, GH is a convex polygon with at most k edges, and can be constructed in O(n log k) time. We then construct the Voronoi diagram of the edges of GH (in O(k log k) time); this allows us to determine, in O(n log k) time, which of the n discs of D lie in GH’s interior. 137 Next, we characterize the curves containing boundary points of the possible 1centers of uncertain discs. We will need this information later in order to construct the possible 1-centers. Lemma 4.7 If D is a set of uncertain discs, then each curve If i lies on the arm of a hyperbola whose foci are of and oi ; and each curve IIij lies on the boundary of a disc. Proof. By definition, every point p ∈ If i is the center of a disc that has outsidetangent Df and inside-tangent Di . It is easy to show that p, of , and oi satisfy the equation of a hyperbola of the stated type. Since the regions are discs, the Minkowski average of any pair of regions is also a disc; hence, each curve IIij lies on a circle. Type III curves for uncertain discs are more complicated. In the following, we show how each such curve can be represented as a system of polynomial equations. Numerical methods can then be applied to this system to generate points on the curve. Let p be a point on curve IIIf ij , Cp be the disc centered at p corresponding to the definition for IIIf ij , and a (resp., b) the critical point where the boundaries of Cp and Di (resp., Dj ) intersect. If we place the axes so that of is at the origin (Figure 4.9), then: d(a, oi ) − ri = 0 (4.1) d(b, oj ) − rj = 0 (4.2) 1 p = (a + b) 2 (4.3) d(p, 0) = d(a, b)/2 + rf . (4.4) Substituting (4.3) into (4.4) leads to 1 1 (ax + bx )2 + (ay + by )2 = 4 4 1 1 (ax − bx )2 + (ay − by )2 + d(a, b) · rf + rf2 ; 4 4 138 and since (u + v)2 − (u − v)2 = 4uv, this can be rearranged to yield a · b − d(a, b) · rf − rf2 = 0 , (4.5) where a · b is the scalar product of a and b. Equations (4.1), (4.2), and (4.5) now form a system of polynomial equations defining pairs of points (a, b) whose midpoints lie on IIIf ab . IIIf ij a p oi b oj Di Dj Cp Df Figure 4.9: Parameterization of IIIf ib . 4.3 Query Algorithm In this section, we present an efficient algorithm to determine whether or not a particular point is a 1-center of a feasible set of a set of uncertain discs. Later, we will extend the algorithm to the case where the uncertain regions are axis-aligned rectangles. To gain some intuition into how the algorithm will work, we will require some additional terminology. See Figure 4.10. Let q be the query point, and Cqr be the disc with radius r centered at q. For each boundary point x of Cqr , let x e denote its antipodal point. Let Df be the farthest region of D from q, r∗ be its distance from 139 q, and pr be the point where the boundary of Cqr crosses the ray from q through the closest point of Df .1 pr Cqr Df pr ∗ q per Figure 4.10: Terminology related to query point q. As an optimization step, we will use the algorithm of Lemma 4.6 to remove any regions from D that do not contribute possible 1-centers. Lemma 4.8 Df must be one of the k regions of D that is not interior to GH. Proof. Let Di be any region of D that is interior to GH, and let L be the line that is →, (ii) left-tangent to D , and (iii) has segment qo strictly (i) perpendicular to ray − qo i i i to its left. By the proof of Theorem 4.5, some region of D lies strictly to the right of L, and is hence farther from q than Di . The intersection and support conditions of smallest bounding discs can be extended to the situation where the points lie within a set of regions of uncertainty. The intersection condition extends directly: every smallest bounding disc of a feasible set of D must intersect every region of D. Lemma 4.9 If r∗ = 0, then q ∈ C◦ ; otherwise, a particular Cqr satisfies the intersection condition of smallest bounding discs if and only if r ≥ r∗ . Proof. If r∗ = 0, then every region of D must intersect q. Our general position assumption then implies that there exists some feasible set of points very close to 1 If r∗ = 0, then this ray is not well-defined; but as we will see later, this will not be a problem. 140 q whose 1-center is at q. If r∗ > 0, then only if r ≥ r∗ will disc Cqr intersect Df (and every region of D). Extending the support condition is more challenging, since a particular uncertain region can contribute at most one point to the critical subset of each smallest bounding disc. To address this issue, we require some additional definitions. If u and v are two boundary points of disc C, then we denote the boundary arc from u ccw to v by arc(u, v). The angle of an arc is the angle subtended by the arc. If the boundary of C intersects at least two regions of D, and s is a boundary point of C, then we define the gap at s to be the smallest arc(u, v) of this boundary that contains s and for which there are distinct regions Di , Dj ∈ D such that u ∈ Di and v ∈ Dj . If the boundary of C intersects fewer than two regions of D, then we will consider the gap at s to be the entire boundary of C, and to have an angle of 2π. A gap is large if its angle is greater than π, and small otherwise. In Figure 4.11, the gap at u is small, while those at v and w are large. u v w Figure 4.11: Examples of gaps at disc boundary points. Lemma 4.10 Point q lies within C◦ iff there exists some Cqr , with r ≥ r∗ , whose every boundary point lies within a small gap. Proof. If q ∈ C◦ , then there exists r ≥ r∗ and a feasible set whose smallest bounding disc is Cqr . Since Cqr satisfies the support condition, every point on its boundary lies within an arc with angle not greater than π whose endpoints lie within distinct regions of D, and hence lies within a small gap. Suppose that for every r ≥ r∗ , at least one point p on the boundary of Cqr lies within a large gap (Lemma 4.9 implies that we can ignore any r < r∗ ). If p lies 141 within some region of D, then since p is within a large gap, no other region of D intersects the boundary of Cqr ; whereas if p is not within any region of D, then p must lie within a boundary arc of Cqr whose angle exceeds π, and which does not intersect any regions of D. We can conclude that no feasible smallest bounding disc of D centered at q can have radius r and satisfy the support condition. Since this holds for all r ≥ r∗ , q ∈ / C◦ . We will now assume that D is a set of uncertain discs. Lemma 4.10 does not immediately motivate an algorithm for determining if q ∈ C◦ , since we cannot afford to examine every point on every disc Cqr to see if it lies within a large gap. Fortunately, we can reduce our work by exploiting a close relationship between query point q, the support condition of feasible smallest bounding discs, and the farthest disc Df from q. Lemma 4.11 Point q lies within C◦ iff there exists some Cqr , with r ≥ r∗ , such that the gaps at pr and per are both small. Proof. We can use the proof of Lemma 4.10 to prove that if q ∈ C◦ , then some Cqr exists (where r ≥ r∗ ) with small gaps at pr and per . To prove the converse, suppose some Cqr exists such that the gaps at pr and per are both small, where r ≥ r∗ . If r∗ = 0, we are done (by Lemma 4.9). Otherwise, we will show that q ∈ C◦ by finding a critical subset for some disc centered at q whose radius is not less than r∗ . Let arc(a, b) be the gap at pr , arc(c, d) be the gap at per , and δ(x) denote the disc of D contributing gap endpoint x ∈ {a, b, c, d}. Since the gaps are minimal, the points (pr , b, c, per , d, a, pr ) must appear in ccw sequence on Cqr ’s boundary. Without loss of generality, there are five cases to consider. Case 1: Each of the four points lie in distinct discs of D. Three of the four points must form a critical subset for Cqr . Case 2: δ(a) = δ(c). Either eb or de must lie in δ(a), which implies that (b, eb) e is a critical subset for C r . or (d, d) q Case 3: δ(a) = δ(d). Either arc(d, a) or pr per must lie within δ(a). In the former case, some (x ∈ arc(d, a), b, c) is a critical subset for Cqr . In the latter 142 case, q ∈ δ(a), which (since r∗ > 0 and q ∈ / Df ) implies that δ(a) 6= Df , and ∗ (pr∗ ∈ Df , per∗ ∈ δ(a)) is a critical subset for Cqr . Case 4: δ(a) = δ(c) and δ(b) = δ(d). This can be dealt with using the same argument as case 2. Case 5: δ(a) = δ(d) and δ(b) = δ(c). If arc(d, a) does not lie within δ(a), or if arc(b, c) does not lie within δ(b), we can apply the argument of case 3. Otherwise, some (x ∈ arc(d, a), y ∈ arc(b, c)) must be a critical subset for Cqr . Lemma 4.11 does not hold for arbitrary antipodal pairs on the boundary of Cqr . For example, in Figure 4.12, neither u nor u e lie within large gaps, yet q is not within C◦ (D). uf q u Figure 4.12: Lemma 4.11 does not hold for u, u e. We now have the tools necessary to construct an efficient algorithm. Our first step is to apply the algorithm of Lemma 4.6 to remove any discs of D that lie within GH. Our second step will be to determine Df , and from it, the value of r∗ . If r∗ is zero, then by Lemma 4.9, we return true. Otherwise, we perform additional steps to determine whether there exists an r ≥ r∗ such that pr and per both lie within small gaps. Per Lemma 4.11, we return true iff such an r is found. We will construct indicator gap functions for both pr and per . The function for pr , for example, is a function that indicates, for each r ≥ r∗ , whether the gap for pr is large. We will actually derive the indicator gap function for a point by combining left and right half gap functions, each representing the angle of the gap −→∗ . lying to the left or right of ray − qp r 143 Our algorithm consists of the following steps: 1. Remove from D all but the k discs that may contribute possible 1-centers. 2. Determine Df and r∗ . If r∗ = 0, return true. 3. Construct left and right half gap functions for pr . 4. Combine half gap functions into indicator gap function for pr . 5. Repeat steps 3 and 4 for per . 6. Use the indicator gap functions to see if pr and per both lie in small gaps for some r ≥ r∗ . If so, return true; else, return false. Let us investigate the third step of the algorithm in more detail. We will consider only the left half gap function for pr , as the others are symmetric. We construct a disc curve for each Di ∈ D as follows. Let B be the portion of Di lying −→∗ . If B = ∅, then we omit D . Otherwise, the disc curve is the to the left of − qp r i → (see Figure 4.13). portion of the boundary of B that lies to the right of − qo i pr ∗ Df q Figure 4.13: Disc curves of left half gap function for pr . →. Consider the homeomorphism H between Let Φ be the polar angle of − qo f polar coordinate points (r ≥ r∗ , Φ ≤ θ ≤ Φ + π) and cartesian coordinate points (r ≥ r∗ , 0 ≤ θ ≤ π). If we apply H to the disc curves, then the left half gap 144 function corresponds to the lower envelope, or 0-level, of the transformed disc curves. By definition, the boundary points of the gap for pr must lie in distinct discs of D. This complicates matters, since it implies that the (full) gap function for pr cannot be defined simply as the sum of the left and right half gap functions. To remedy this situation, we construct both the 0-level and 1-level of the arrangement of the transformed disc curves. Note that for each value of r, a particular disc Di will contribute a value to at most one of these levels. Hence, if we have these levels for both the left and right half gap functions, and we have recorded which disc has induced each edge of these levels, we can determine the value of the full gap function. To clarify, for the left half gap function, we define three functions: θ0L (r), which returns the value of θ corresponding to the point of the 0-level at distance r from q; θ1L (r), which returns the corresponding value for the 1-level; and δ L (r), which returns the disc of D corresponding to θ0L (r). We define analogous functions θ0R (r), θ1R (r), and δ R (r) for the right half gap function. The value of the full gap function for a particular r ≥ r∗ is now θ0L (r) + θ0R (r) (if δ L (r) 6= δ R (r)) or min(θ0L (r) + θ1R (r), θ1L (r) + θ0R (r)) (otherwise). The fourth step of the algorithm, combining the left and right half gap functions for pr into a full gap function, can be divided into two stages: a merge stage, and a convert stage. For the merge stage, we first take the four levels generated in the previous step (the 0-level and 1-level of the left and right half gap functions) and combine them to form the three functions θ0L (r) + θ0R (r), θ0L (r) + θ1R (r), and θ1L (r) + θ0R (r). We modify the first of these so that for edges where the discs associated with the left and right 0-levels are the same (and are hence illegal), we return 2π. We next merge these three functions into the full gap function: for a particular r, it returns the gap angle for pr as the minimum of the three constituent functions’ values. In the convert stage, we convert the current half gap function (a continuous function of r) into a function that indicates whether the gap for pr is large or small. The final step of the algorithm is to perform a simultaneous scan of the indicator gap function sequences for pr and per , to see if there exists a value r where neither gap angle exceeds π. If one is found, we can return true immediately; otherwise, 145 when the scan is complete, we return false. Let us determine the running time of our algorithm. By Lemma 4.6, step 1 requires O(n log k) time, where k is the number of discs not interior to GH. Step 2 clearly runs in O(n) time. For step 3, we require the following lemma. Lemma 4.12 The 0-level of a half gap function of k uncertain discs has O(k) complexity. Proof. Let D be a set of k discs. Without loss of generality, we assume we are dealing with the left half gap function of pr , and that of lies on the positive x-axis. In addition, for ease of exposition, we will assume that no disc (other than Df ) intersects this axis. Let S be a sequence of representative interior points of each edge of the 0-level, ordered by their distance from the origin. Assume, by way of contradiction, that there exists a subsequence (a1 , b2 , a3 , b4 ) of S, where a1 and a3 are points on the curve for a ∈ D, and b2 and b4 are points on the curve for some other disc b ∈ D. By proving that no such subsequence can exist, we will show that S is a Davenport-Schinzel sequence of length λ2 (k), which implies that |S| = O(k). bb44 a3 a1 b2 R q r1 r2 r3 r4 Figure 4.14: Lemma 4.12. For each point of the subsequence, we construct an arc centered at q that extends from the positive x-axis to the point. Let the radii of these four arcs be r1 < r2 < r3 < r4 (Figure 4.14). Consider the region R that is bounded by the 146 x-axis, segment a1 a3 , and the arcs for r1 and r3 . Note that b2 must lie within R, otherwise it would not lie on a 0-level edge. It follows that segments a1 a3 and b2 b4 intersect. Now, these two segments are chords of a and b respectively, and it can be shown that if two such chords intersect, then at least one of the chord’s endpoints must lie within the other chord’s disc. This is a contradiction, since the endpoint in question cannot lie in the interior of a 0-level edge. Lemma 4.13 The 1-level of a half gap function of k uncertain discs has O(kα(k)) complexity. Proof. Observe that we can add almost-vertical rays extending from the endpoints of each transformed disc curve to produce an unbounded curve that yields the same 0- and 1-levels as the original transformed curves. Note also that each pair of such modified curves will intersect in at most two points. Hence, these curves can be viewed as pseudoparabolas [69], which implies that the 1-level of the half gap function has O(kα(k)) complexity [3]. Lemma 4.12 implies that we can construct the 0-level for each half gap function in O(k log k) time, by modifying the standard divide-and-conquer algorithm of [6]. If we then clip away the edges forming the 0-level from the transformed disc curves, then the lower evelope of the curves that remain form the 1-level of the original transformed disc curves. This clipping can be performed in linear time, and Lemma 4.13 implies that we can use Hershberger’s algorithm [31] to construct the 1-level in O(k log k) time. The merge stage of step 4 can be performed in time linear in the complexity of the 0- and 1-levels, and hence runs in O(kα(k)) time. The convert stage running time is linear in the size of its output, which is also O(kα(k)): Lemma 4.14 Every indicator gap function of k uncertain discs has O(kα(k)) complexity. Proof. The (output) indicator gap function is derived from the (input) original gap function. Since the input is represented by a sequence of length O(kα(k)), we prove that the output sequence has the same size by showing that each input element generates O(1) output elements. 147 We will represent the input gap function by a sequence of triples ((r1 , δ1L , δ1R ), (r2 , δ2L , δ2R ), . . .), where δiL (resp., δiR ) is the disc inducing the left (resp., right) half gap angles for values ri ≤ r < ri+1 . Our task is to convert this sequence into an indicator gap function, a sequence of pairs ((r1 , f1 ) . . .) where each fi indicates whether the gap is large for values immediately following ri . We can perform this task by iterating over the original gap function. Consider a particular step in this iteration, where the current triple is (ru , Di , Dj ), and the following triple is (rv , ·, ·). At this point, we must determine those values of r between ru and rv where the gap function equals π. Geometrically, this is equivalent to finding points on the boundaries of discs Di and Dj whose midpoint is q. Each such r corresponds to a point where the gap angle stops or starts being greater than π. Di a q Dj ru r rv b Figure 4.15: Lemma 4.14. Figure 4.15 shows one triple being processed by the convert stage. Every value 148 r that generates a new output element must satisfy the following equations: d(a, oi ) = ri d(b, oj ) = rj bx = −ax by = −ay r = d(q, a) These equations can be simplified to yield ay = c1 ax + c2 , the equation of a line with constant coefficients c0 = oi .y + oj .y c1 = −(oi .x + oj .x)/c0 c2 = (oi .x2 + oi .y 2 − ri2 − oj .x2 − oj .y 2 + rj2 )/c0 . Since a must be one of the two possible points of intersection of this line with Di , and each r corresponds to a unique point a, there are at most two values of r (and hence at most two output elements) generated per triple. The proof of Lemma 4.14 implies that step 4 runs in O(kα(k)) time. Step 5 can be charged to steps 3 and 4. Step 6 involves a simple simultaneous traversal of the two indicator gap function sequences, and hence has the same running time as step 4. The running time of the complete algorithm is thus dominated by step 1, and we can make the following claim. Theorem 4.15 If D is a set of n uncertain discs, then it can be determined if a point q is within C◦ (D) in O(n log k) time, where k is the number of discs of D that do not lie in the interior of GH(D). 149 4.3.1 Uncertain Axis-Aligned Rectangles In this section, we show how the query algorithm can be adapted to the case where the uncertain regions are axis-aligned rectangles. We will require the following analog to Lemma 4.11 (its proof can be found in Appendix A). Lemma 4.16 If D is a set of uncertain axis-aligned rectangles, then point q lies within C◦ iff there exists some Cqr , with r ≥ r∗ , such that the gaps at pr and per are both small. The structure of the query algorithm is the same as before, but some of the steps are different. Let us consider step 3: constructing the left and right half gap functions for pr . As we did for the case of uncertain discs, we will consider only the left half gap function for pr . For each rectangle Di , we construct a set of boundary edges (analogous to the disc curves for uncertain discs). These are −→∗ are segments from the boundary of Di , after portions lying to the right of − qp r removed; see Figure 4.16. Df pr ∗ q Figure 4.16: Boundary edges for axis-aligned rectangles. Lemma 4.17 The 0- and 1-levels of a half gap function of k uncertain axis-aligned rectangles have O(kα(k)) complexity. Proof. There are O(k) boundary edges induced by a set of k such rectangles, and since each such edge is a line segment, both levels of their arrangement have O(kα(k)) complexity [3]. 150 The 0- and 1-levels are derived from an arrangement of line segments. This implies that we can construct the 0-level in O(k log k) time, by using Hershberger’s algorithm [31]. To construct the 1-level, we take the same approach as before: we clip away the portions of the line segments comprising the 0-level (in linear time), and construct the lower envelope of the remaining segments. By Lemma 4.17, we can use Hershberger’s algorithm to construct the 1-level in O(k log k) time. To prove that the algorithm for uncertain axis-aligned rectangles has the same running time as for uncertain discs, all that remains is to derive an analog to Lemma 4.14: Lemma 4.18 Every indicator gap function of k uncertain axis-aligned rectangles has O(kα(k)) complexity. Proof. By Lemma 4.17, the gap function for such rectangles has O(kα(k)) complexity. Hence, it will suffice to show that each element of this gap function induces O(1) elements in the corresponding indicator gap function. We can assume that the input gap function is represented by a sequence of triples ((r1 , δ1L , δ1R ), (r2 , δ2L , δ2R ), . . .), where δiL (resp., δiR ) is the boundary segment inducing the left (resp., right) half gap angles for values ri ≤ r < ri+1 . Suppose the gap function for some r is equal to π, where (ri , δiL , δiR ) is an element in this sequence, and ri ≤ r < ri+1 . Observe that one of the boundary segments will be vertical, and the other will be horizontal.2 Without loss of generality, we will assume that δrL is horizontal, as in Figure 4.17. Then we have: ay = c1 bx = c2 bx = −ax by = −ay r = d(q, a) , 2 Except when the two boundary segments lie on parallel lines equidistant from q; but in that case, the gap function for every r in the interval is equal to π, including the endpoints, and so does not contribute to the indicator gap function sequence. 151 where c1 and c2 are constants. These equations simplify to r= q c21 + c22 ; hence at most one value of r is generated for each element of the input gap function sequence. a q r b Figure 4.17: Lemma 4.18. Since we have shown that the running time of each step of the algorithm is the same as it was in the case for uncertain discs, we can make the following claim. Theorem 4.19 If D is a set of n uncertain axis-aligned rectangles, then it can be determined if a point q is within C◦ (D) in O(n log k) time, where k is the number of rectangles of D that do not lie in the interior of GH(D). 4.4 Construction Algorithm In this section, we present our algorithm for constructing the possible 1-centers of a set of uncertain discs. By Theorem 4.4, the boundary of C◦ lies on type I, II, and III curves. Our algorithm uses a plane sweep to construct the arrangement of these curves, then extracts the subset of the edges of the arrangement that form the boundary of C◦ . As we did for the previous algorithm, we first remove from D all but the k discs that do not lie in the interior of GH. For efficiency, we would like to minimize the number of possible intersections between the curves involved. Since 152 each curve IIIf ij lies in Df ’s cell within the farthest disc Voronoi diagram of D, if we restrict our attention to the portion of C◦ lying within the cell for Df , we need only consider intersections between pairs of k2 possible type III curves (i.e., when f is fixed), as opposed to intersections between pairs of k3 type III curves (i.e., when f can vary) within the complete C◦ . Since the number of intersections is a function of the square of the number of curves involved, this will improve our guaranteed worst-case running time by a factor of k, as we shall see. Our algorithm has the following steps: 1. Remove from D all but the k discs that may contribute possible 1-centers. 2. Construct V , the farthest disc Voronoi diagram of D. 3. For each cell Rf ∈ V : (a) Generate every type I, II, and III curve that may contain boundary points of C◦ within Rf . (b) Construct A, the arrangement of these curves and the boundary of Rf . (c) For each edge E of A, choose point u in the interior of E, and points v and w lying very close to and on opposite sides of E. If u is within Rf , and exactly one of v and w lies within C◦ , output E. Theorem 4.4 implies that the interior points of any particular cell of A are either all within, or all outside of, C◦ . Hence, we need only to determine which edges of A are boundary edges of C◦ . As we process each edge E of each cell Rf , we use the test of point u to verify that E lies within the cell. If this test is passed, then we verify that E is a boundary of C◦ by testing points v and w. We can conclude that the algorithm constructs the boundary edges of C◦ . The running time of step 1 is O(n log k) (Lemma 4.6). Step 2 takes O(k log k) time [60]. For a particular cell Rf ∈ V , there are O(k) candidate hyperbolic arcs containing type I curves of the form If ∗ that may intersect the cell, and k2 pairs of discs (Di , Dj ) associated with curves IIij and IIIf ij that may intersect the cell. The total number of curves generated in step 3(a) for all O(k) cells of V is thus O(k 2 ). 153 Since each curve is a polynomial curve of bounded degree, the number of intersections occurring between pairs of curves is at most O(k 4 ). If we assume that (i) each curve can be split into O(1) pieces, each monotonic in x, in O(1) time, and (ii) all intersections between any pair of these pieces can be determined in O(1) time, then A can be constructed in step 3(b) by using the plane sweep algorithm of Clarkson and Shor [16] or Mulmuley [53] in expected O(k 4 ) time. We can determine which of points u, v, and w lie within C◦ (D), in O(k log k) time, by using the algorithm of Section 4.3. Since there are O(k 4 ) edges in A, step 3(c) takes O(k 5 log k) time; and since there are O(k) cells in V , the total expected time spent in step 3 is O(k 6 log k). The running time of the algorithm is thus dominated by steps 1 and 3, and we can make the following claim. Theorem 4.20 The possible 1-centers of n uncertain discs D can be constructed in expected O((n + k 6 ) log k) time, where k is the number of discs of D that do not lie in the interior of GH(D). 4.4.1 Uncertain Axis-Aligned Rectangles We now consider the case where each Di ∈ D is an (axis-aligned) rectangle (Figure 4.18). Observe that since each rectangle has nonzero area, a rectangle can be outside-tangent to a disc at either a vertex or a point along an edge, whereas it can be locally inside-tangent to a disc only at a vertex. Note also that with these rectangles, local inside tangency does not imply ‘normal’ inside tangency. Lemma 4.21 If D is a set of axis-aligned rectangles, then each curve If i lies on O(1) lines and parabolas, and each curve IIij lies on the boundary of a rectangle. Proof. By definition, each point s ∈ If i is the center of a disc that has outsidetangent (at either a corner or an edge) Df and locally inside-tangent (at a corner) Di . Each corner/corner combination lies on the bisectors of the corners, whereas each edge/corner combination lies on a parabola induced by the edge and corner. The claim for IIij follows from the fact that the Minkowski average of a pair of axis-aligned rectangles is also a rectangle. 154 Figure 4.18: Possible 1-centers of uncertain rectangles. Lemma 4.22 If D is a set of axis-aligned rectangles, then every point on the boundary of C◦ that is not on a type I or type II curve will lie on one of a small number of lines or parabolas. Proof. Let p be a point on the boundary of C◦ that is not on a type I or type II curve. By Theorem 4.4, there must exist a feasible set S of D with a critical subset of boundary sites {si , sj }, and B(S) must have outside-tangent some Df . Without loss of generality, we can assume that si .x ≥ sj .x. Observe that si and sj cannot both lie on the left sides of their respective rectangles, otherwise p would lie on Iij . If si lies in the interior of its rectangle’s left side, and sj lies in the interior of its rectangle’s right side, then it is possible to move si and sj in opposite directions along these sides, to form a new feasible set S 0 with the same 1-center as S, but with a larger (smallest) bounding disc (Figure 4.19). It is then possible to move s0i and s0j directly away from each other to yield a feasible set S 00 whose 1-center remains at p but that has no boundary sites, contradicting Lemma 4.2. Symmetric arguments hold for the cases where si and sj both lie on the right, top, or bottom sides of their respective rectangles, or when si and sj both lie in the interior of horizontal sides. Hence it must be the case that si and sj lie on sides of their respective rectangles that have opposite orientations (i.e., one vertical, and one horizontal).3 3 This includes the rectangle corners, which are intersections of edges of both orientations. 155 s00i Di s0i si p sj s0j Dj s00j Df Figure 4.19: Lemma 4.22. Df can be outside-tangent B(S) in two ways. First, the point of tangency can be at a corner (xf , yf ) of Df . In this case, we can assume without loss of generality that (i) si lies on a horizontal side of Di at position (x, yi ), (ii) sj lies on a vertical side of Dj at position (xj , y), where xf , yf , xj , and yi are constants (Figure 4.20). Then 1 px = (x + xj ) 2 1 py = (yi + y) 2q q 1 2 2 (px − xf ) + (py − yf ) = (x − xj )2 + (yi − y)2 . 2 This can be simplified to yield py = c1 px + c2 , 156 the equation of a line with constant coefficients c1 = (xj − xf )/(yf − yi ) c2 = (x2f + yf2 − x2j + 1)/(2(yf − yi )) . Dj Df sj = (xj , y) (xf , yf ) B(S) C ◦(D) p IIIf ij Di si = (x, yi) Figure 4.20: IIIf ij is a line. Alternatively, Df can be outside-tangent B(S) along a side of Df . Then without loss of generality, we can make the same assumptions as in the previous case, except that now B(S) is outside-tangent a horizontal side of Df lying on the line y = yf (Figure 4.21). Our system of equations now becomes 1 px = (x + xj ) 2 1 py = (yi + y) 2q 1 |py − yf | = (x − xj )2 + (yi − y)2 . 2 157 This can be simplified to yield py = c1 p2x + c2 px + c3 , the equation of a parabola with constant coefficients c1 = 1/(2(yi − yf )) c2 = −2xj c1 c3 = (x2j − yf2 − 1)c1 . Df Dj (px , yf ) B(S) sj = (xj , y) C ◦ (D) IIIf ij p Di si = (x, yi ) Figure 4.21: IIIf ij is a parabola. As there is a small number of possible choices of sides of Di and Dj and sides or corners of Df associated with IIIf ij , every such boundary point p must lie on some small number of lines and parabolas. The preceding results allow us to make the following claim. Theorem 4.23 The possible 1-centers of n uncertain axis-aligned rectangles D can be constructed in expected O((n + k 6 ) log k) time, where k is the number of rectangles of D that do not lie in the interior of GH(D). 158 4.5 Closing Remarks In this chapter, we studied a variant of a classical facility location problem, the smallest bounding disc of a set of points, in which each site lies in a region of uncertainty. We introduced a new geometric figure: the possible 1-centers of uncertain regions. In addition to presenting two algorithms related to this figure, we also demonstrated a connection between this figure and the guaranteed hull of the regions. We have identified a number of open problems concerning possible 1-centers of uncertain discs. We do not, as yet, have bound on the complexity of these structures that is tighter than the O(n5 ) bound implied by our construction algorithm (Section 4.4). We also have not been able to prove that C◦ is simply connected (i.e., has no holes), although we conjecture that this is true. Our construction algorithm can be simplified if this is the case. Figure 4.22: Possible bounding disc, guaranteed bounding disc (uncertain discs). Figure 4.23: Possible bounding disc, guaranteed bounding disc (uncertain axis-aligned rectangles). 159 Two additional problems related to possible 1-centers that are candidates for future research involve what we call the possible bounding disc and the guaranteed bounding disc of uncertain regions (Figures 4.22 and 4.23). We define the first as the union of all feasible smallest bounding discs (in contrast to the origins of these discs, which form the regions’ possible 1-centers). We define the second as the intersection of all feasible smallest bounding discs. An applet exploring some of the geometric objects described in this chapter can be found at http://www.cs.ubc.ca/∼jpsember/ocent.html. The applet (Figure 4.24) includes the following sample files: 1. p1c.dat: Generates possible 1-centers of a set of uncertain discs. The program uses a simple (and slow) random sampling process to do this, instead of the algorithm described in Section 4.4. 2. p1c rects.dat: Generates possible 1-centers of a set of uncertain axisaligned rectangles. 3. psbd.dat: Constructs the possible bounding disc of a set of uncertain discs or axis-aligned rectangles. 4. g1c.dat: Constructs the guaranteed bounding disc for a set of uncertain discs or axis-aligned rectangles. Figure 4.24: Applet for Chapter 4. 160 Chapter 5 Conclusion and Future Research In this thesis, we considered problems from computational geometry involving point sets, and studied variants of the problems where each point’s location was replaced by a region of uncertainty. We can often view the output of each of the original problems as being a partition of the plane into two sets: points that satisfy some property (e.g., lie within the convex hull of the points), and those that do not.1 The partition of the plane associated with the output of our modified problems includes a third set: points that we cannot make any claims about with respect to the property. We have referred to the points within the first of these sets as a guaranteed object, and to the points in the first or third of these sets as a possible object. In the introduction, we showed that imprecision is an unavoidable factor in most geometric algorithms, whether as a result of imprecision in the input, or as imprecision introduced by the representation and manipulation of numbers by the algorithm. By classifying the output of an algorithm as either a guaranteed or possible object, an algorithm may more precisely represent those properties that can be safely inferred from the input data. In each of the next three chapters, we applied this technique to a different fundamental computational geometry problem. In Chapter 2, we focused on Voronoi diagrams of uncertain sites. We identi1 Unless the output of the algorithm is combinatoric in nature, in which case the set being partitioned is also combinatoric. For example, the Delaunay triangulation of a point set is a subset of the edges of the complete graph of the sites. 161 fied as an area of future research the possibility of generating order-k guaranteed Voronoi diagrams, and discussed some of the issues likely to be encountered while doing so. We also left as future work the problem of constructing guaranteed Voronoi diagrams for non-disjoint uncertain polygons. We defined guaranteed variants of two combinatoric structures related to Voronoi diagrams: guaranteed Delaunay edges, and guaranteed Gabriel edges. We briefly discussed minimum spanning trees and relative neighborhood graphs, two other combinatoric structures associated with point sites that may yield interesting variants when the sites are imprecise. In Chapter 3, we addressed two variants of the convex hull problem involving uncertain sites: guaranteed hulls, and possible hulls. We highlighted the fact that existing algorithms for generating guaranteed hulls are closely tied to hull bitangents. In some situations this may be undesirable, since there may be a large number of hull bitangents that do not contribute to the guaranteed hull. Finding a guaranteed hull algorithm with a reduced dependency on these hull bitangents is thus an open problem. We also left as an open problem the task of tightening the O(nα(n)) bound on the number of hull bitangents of uncertain polygons. In Chapter 4, we considered problems related to smallest bounding discs of imprecise points. Some fundamental questions regarding possible 1-centers of uncertain discs (or uncertain axis-aligned rectangles) remain: What is the complexity of this structure? And, is it simply connected? In each of the three preceding chapters, we have referred to software programs that demonstrate the research presented in this thesis. The programs are written in Java, and are built upon a geometry editor and algorithm debugging framework that can be found at http://www.cs.ubc.ca/∼jpsember/testbed.html. This framework lets users create, manipulate, and experiment with sets of geometric objects in a visual manner, to gain insight into a wide range of computational geometry problems. Many of the results we presented in this thesis were sparked by observations made while using this software. 162 Bibliography [1] M. Abellanas, F. Hurtado, and P. A. Ramos. Tolerance of geometric structures. In Proceedings of the 6th Canadian Conference on Computational Geometry, pages 250–255, 1994. → pages 16 [2] M. Abellanas, F. Hurtado, and P. A. Ramos. Structural tolerance and Delaunay triangulation. Information Processing Letters, 71(5-6):221–227, 1999. → pages 16 [3] P. K. Agarwal, B. Aronov, T. M. Chan, and M. Sharir. On levels in arrangements of lines, segments, planes, and triangles. Discrete and Computational Geometry, 19(3):315–331, 1998. → pages 147, 150 [4] M. E. Ali, E. Tanin, R. Zhang, and K. Ramamohanarao. Probabilistic Voronoi diagrams for processing moving nearest neighbor queries. Technical report, University of Melbourne, Melbourne, Australia, 2009. → pages 15 [5] T. Asano, J. Matoušek, and T. Tokuyama. Zone diagrams: Existence, uniqueness, and algorithmic challenge. In Proceedings of the 18th Annual ACM-SIAM Symposium on Discrete Algorithms, volume 37, pages 756–765, 2007. → pages 15 [6] M. Atallah and C. Bajaj. Efficient algorithms for common transversals. Information Processing Letters, 25(2):87–91, 1987. → pages 85, 147 [7] F. Aurenhammer. Voronoi diagrams—a survey of a fundamental geometric data structure. ACM Computing Surveys, 23(3):345–405, 1991. → pages 12 [8] F. Aurenhammer, G. Stöckl, and E. Welzl. The post office problem for fuzzy point sets. In Proceedings of the International Workshop on Computational Geometry - Methods, Algorithms and Applications, pages 1–11, London, UK, 1991. Springer-Verlag. → pages 14 163 [9] C. B. Barber, D. P. Dobkin, and H. Huhdanpaa. The Quickhull algorithm for convex hulls. ACM Transactions on Mathematical Software, 22(4):469–483, 1996. → pages 71 [10] B. Bhattacharya, S. Jadhav, A. Mukhopadhyay, and J. Robert. Optimal algorithms for some intersection radius problems. Computing, 52(3): 269–279, 1994. → pages 129 [11] J. Boisson, A. Cérézo, O. Devillers, J. Duquesne, and M. Yvinec. An algorithm for constructing the convex hull of a set of spheres in dimension d. Computational Geometry: Theory and Applications, 6(2):123–130, 1996. → pages 99, 100, 125 [12] K. Q. Brown. Voronoi diagrams from convex hulls. Information Processing Letters, 9(5):223–228, 1979. → pages 11 [13] K. Buchin, M. Löffler, P. Morin, and W. Mulzer. Delaunay triangulation of imprecise points simplified and extended. In Proceedings of the 11th International Symposium on Algorithms and Data Structures, pages 131–143, Berlin, Heidelberg, 2009. Springer-Verlag. → pages 16 [14] E. W. Chambers, A. Erickson, S. P. Fekete, J. Lenchner, J. Sember, V. Srinivasan, U. Stege, S. Stolpner, C. Weibel, and S. Whitesides. Connectivity graphs of uncertainty regions. In Proceedings of the 21st International Symposium on Algorithms and Computation, 2010. → pages 68 [15] T. M. Chan. Optimal output-sensitive convex hull algorithms in two and three dimensions. Discrete and Computational Geometry, 16:361–368, 1996. → pages 72, 75, 86, 93, 94 [16] K. L. Clarkson and P. W. Shor. Applications of random sampling in computational geometry, ii. Discrete and Computational Geometry, 4: 387–421, 1995. → pages 154 [17] H. Davenport and A. Schinzel. A combinatorial problem connected with differential equations. American Journal of Mathematics, 87(3):684–694, 1965. → pages 6 [18] M. de Berg, M. van Kreveld, M. Overmars, and O. Schwarzkopf. Computational Geometry: Algorithms and Applications. Springer-Verlag, second edition, 2000. → pages 31 164 [19] O. Devillers. Delaunay triangulation of imprecise points, preprocess and actually get a fast query time. In XIV Spanish Meeting on Computational Geometry, 2011. → pages 16 [20] R. L. Drysdale and D. T. Lee. Generalized Voronoi diagrams in the plane. In Proceedings of the 16th Annual Allerton Conference on Communications, Control and Computing, pages 833–842, 1978. → pages 14 [21] A. Edalat and A. Lieutier. Foundation of a computable solid modelling. Theoretical Computer Science, 284(2):319–345, 2002. → pages 73 [22] A. Edalat, A. Lieutier, and E. Kashefi. The convex hull in a new model of computation. In Proceedings of the 13th Canadian Conference on Computational Geometry, pages 93–96, 2001. → pages 73 [23] A. Edalat, A. A. Khanban, and A. Lieutier. Computability in computational geometry. In Conference on Computability in Europe, pages 117–127, 2005. → pages 73, 74 [24] J. S. Ely and A. P. Leclerc. Correct Delaunay triangulation in the presence of inexact inputs and arithmetic. Reliable Computing, 6(1):23–38, 2000. → pages 15 [25] E. Ezra and W. Mulzer. Convex hull of imprecise points in o(n log n) time after preprocessing. In Proceedings of the 27th European Workshop on Computational Geometry, 2011. → pages 74 [26] S. Fortune. A sweepline algorithm for Voronoi diagrams. In Proceedings of the 2nd Annual Symposium on Computational Geometry, pages 313–322, New York, NY, USA, 1986. ACM. → pages 11, 14, 25, 36, 42, 51 [27] K. R. Gabriel and R. R. Sokal. A new statistical approach to geographic variation analysis. Systematic Biology, 18(3):259–278, 1969. → pages 14 [28] B. Gärtner. Fast and robust smallest enclosing balls. In Proceedings of the 7th Annual European Symposium on Algorithms, pages 325–338, 1999. → pages 127 [29] R. L. Graham. An efficient algorithm for determining the convex hull of a finite planar set. Information Processing Letters, 1(4):132 – 133, 1972. → pages 71, 93 [30] L. Guibas, D. Salesin, and J. Stolfi. Epsilon geometry: Building robust algorithms from imprecise computations. In Proceedings of the Fifth Annual Symposium on Computational Geometry, 1989. → pages 4 165 [31] J. Hershberger. Finding the upper envelope of n line segments in o(n log n) time. Information Processing Letters, 33(4):169–174, 1989. → pages 147, 151 [32] S. Jadhav, A. Mukhopadhyay, and B. Bhattacharya. An optimal algorithm for the intersection radius of a set of convex polygons. Journal of Algorithms, 20(2):244 – 267, 1996. → pages 129 [33] R. A. Jarvis. On the identification of the convex hull of a finite set of points in the plane. Information Processing Letters, 2(1):18–21, 1973. → pages 71 [34] V. Karamcheti, C. Li, I. Pechtchanski, and C. Yap. A Core library for robust numeric and geometric computation. In ACM Symposium on Computational Geometry, Applied Track, 1999. → pages 3 [35] M. I. Karavelas and M. Yvinec. Dynamic additively weighted Voronoi diagrams in 2d. In Proceedings of the 10th Annual European Symposium on Algorithms, pages 586–598, 2002. → pages 41 [36] A. Khanban. Basic Algorithms of Computational Geometry with Imprecise Input. PhD thesis, Imperial College of London, 2005. → pages 15, 16 [37] A. Khanban and A. Edalat. Computing Delaunay triangulation with imprecise input data. In Proceedings of the 15th Canadian Conference on Computational Geometry, pages 94–97, 2003. → pages 15 [38] D. G. Kirkpatrick. Efficient computation of continuous skeletons. In Proceedings of the 20th Annual Symposium on Foundations of Computer Science, pages 18–27, 1979. → pages 12, 14, 36 [39] D. G. Kirkpatrick and J. D. Radke. A framework for computational morphology. In G. T. Toussaint, editor, Computational Geometry, pages 217–248, Amsterdam, Netherlands, 1985. North-Holland. → pages 68 [40] D. G. Kirkpatrick and R. Seidel. The ultimate planar convex hull algorithm. SIAM Journal on Computing, 15(1):287–299, 1986. → pages 72 [41] R. Klein, K. Mehlhorn, and S. Meiser. Randomized incremental construction of abstract Voronoi diagrams. Computational Geometry: Theory and Applications, 3(3):157–184, 1993. → pages 41 [42] D. E. Knuth. The Art of Computer Programming, Volume III: Sorting and Searching. Addison-Wesley, 1973. → pages 11 166 [43] D. T. Lee and B. J. Schachter. Two algorithms for constructing a Delaunay triangulation. In International Symposium on Computer and Information Sciences, 1980. → pages 11 [44] E. H. Lockwood. A Book of Curves. Cambridge University Press, Cambridge, England, 1967. → pages 57 [45] M. Löffler. Data Imprecision in Computational Geometry. PhD thesis, Utrecht University, 2009. → pages 1 [46] M. Löffler and J. Snoeyink. Delaunay triangulation of imprecise points in linear time after preprocessing. Computational Geometry: Theory and Applications, 43(3):234–242, 2010. → pages 16 [47] M. Löffler and M. J. van Kreveld. Largest and smallest convex hulls for imprecise points. Algorithmica, 56(2):235–269, 2010. → pages 7, 74 [48] M. Löffler and M. J. van Kreveld. Largest bounding box, smallest diameter, and related problems on imprecise points. In Proceedings of the 10th Workshop on Algorithms and Data Structures, pages 447–458, 2007. → pages 129 [49] N. Megiddo. Linear-time algorithms for linear programming in R3 and related problems. SIAM Journal on Computing, 12(4):759–776, 1983. → pages 129 [50] A. A. Melkman. On-line construction of the convex hull of a simple polyline. Information Processing Letters, 25(1):11–12, 1987. → pages 109, 122 [51] R. E. Miles. On the homogeneous planar Poisson point process. Mathematical Biosciences, 6:85 – 127, 1970. → pages 14 [52] R. E. Moore. Interval Analysis. Prentice-Hall, 1966. → pages 4 [53] K. Mulmuley. A fast planar partition algorithm, ii. Journal of the ACM, 38 (1):74–103, 1991. → pages 154 [54] T. Nagai and N. Tokura. Tight error bounds of geometric problems on convex objects with imprecise coordinates. In Revised Papers from the Japanese Conference on Discrete and Computational Geometry, pages 252–263, London, UK, 2001. Springer-Verlag. → pages 74, 75, 82, 83, 86, 93 167 [55] T. Nagai, Y. Seigo, and N. Tokura. Convex hull problem with imprecise input. In Revised Papers from the Japanese Conference on Discrete and Computational Geometry, pages 207–219, London, UK, 2000. Springer-Verlag. → pages 74, 76, 83, 85, 93, 100, 101, 125 [56] F. Nielsen and M. Yvinec. An output-sensitive convex hull algorithm for planar objects. International Journal of Computational Geometry and Applications, 8:39–65, 1995. → pages 93 [57] F. P. Preparata. An optimal real-time algorithm for planar convex hulls. Communications of the ACM, 22(7):402–405, 1979. → pages 87 [58] F. P. Preparata and S. J. Hong. Convex hulls of finite sets of points in two and three dimensions. Communications of the ACM, 20(2):87–93, 1977. → pages 71 [59] F. P. Preparata and D. E. Muller. Finding the intersection of n half-spaces in time o(n log n). Theoretical Computer Science, 8:45–55, 1979. → pages 93, 99 [60] D. Rappaport. Computing the furthest site Voronoi diagram for a set of discs (preliminary report). In Proceedings of the Workshop on Algorithms and Data Structures, pages 57–66, 1989. → pages 153 [61] D. Rappaport. A convex hull algorithm for discs, and applications. Computational Geometry: Theory and Applications, 1(3):171–187, 1992. → pages 73, 74, 99 [62] H. Rosenberger. Order-k Voronoi diagrams of sites with additive weights in the plane. Algorithmica, 6(4):490–521, 1991. → pages 14, 49, 51 [63] J. Sember and W. Evans. Guaranteed Voronoi diagrams of uncertain sites. In Proceedings of the 20th Canadian Conference on Computational Geometry, pages 207–210, 2008. → pages 17 [64] M. I. Shamos. Computational Geometry. PhD thesis, New Haven, CT, USA, 1978. → pages 87 [65] M. I. Shamos and D. Hoey. Closest-point problems. In Proceedings of the 16th Annual Symposium on Foundations of Computer Science, pages 151–162, Washington, DC, USA, 1975. IEEE Computer Society. → pages 11, 14, 26, 30, 129 168 [66] M. Sharir. Intersection and closest-pair problems for a set of planar discs. SIAM Journal on Computing, 14(2):448–468, 1985. → pages 12, 14, 24, 43 [67] M. Sharir and P. K. Agarwal. Davenport-Schinzel Sequences and their Geometric Applications. Cambridge University Press, New York, NY, USA, 2010. → pages 6, 83, 84 [68] J. J. Sylvester. A question in the geometry of situation. 1, 1857. → pages 129 [69] H. Tamaki and T. Tokuyama. How to cut pseudo-parabolas into segments. In Proceedings of the Eleventh Annual Symposium on Computational Geometry, pages 230–237, New York, NY, USA, 1995. ACM. → pages 147 [70] T. Tao. Analysis, Volume II. Hindustan Book Agency, 2006. → pages 130 [71] G. T. Toussaint. The relative neighbourhood graph of a finite planar set. Pattern Recognition, 12(4):261 – 268, 1980. → pages 67 [72] G. T. Toussaint. Solving geometric problems with the rotating calipers. In Proceedings IEEE MELECON, 1983. → pages 71, 122 [73] R. Wein. Efficient implementation of red-black trees with split and catenate operations. Technical report, Tel-Aviv University, Tel-Aviv, Israel, 2005. → pages 92 [74] E. Welzl. Smallest enclosing disks (balls and ellipsoids). New Results and New Trends in Computer Science, pages 359–370, 1991. → pages 129 [75] Y. Yang, M. Lin, J. Xu, and Y. Xie. Minimum spanning tree with neighborhoods. In Proceedings of the 3rd International Conference on Algorithmic Aspects in Information and Management, pages 306–316, Berlin, Heidelberg, 2007. Springer-Verlag. → pages 68 [76] C. K. Yap. An o(n log n) algorithm for the Voronoi diagram of a set of simple curve segments. Discrete and Computational Geometry, pages 365–393, 1987. → pages 36 [77] C. K. Yap. Towards exact geometric computation. Computational Geometry, 7(1-2):3 – 23, 1997. → pages 4 169 Appendix A Proof of Lemma 4.16 We will prove that if D is a set of uncertain axis-aligned rectangles, then point q lies within C◦ iff there exists some Cqr , with r ≥ r∗ , such that the gaps at pr and per are both small. As was the case for Lemma 4.11, we can use the proof of Lemma 4.10 to prove that if q ∈ C◦ , then some Cqr exists (where r ≥ r∗ ) with small gaps at pr and per . We prove the converse as follows. Without loss of generality, we will assume that pr lies within quadrant 1. Our proof will manipulate a disc Cq that is centered at q and intersects Df (and hence satisfies the intersection condition of smallest bounding discs). It will then show that either Cq , or some smaller disc centered at q, satisfies the support condition as well. This will allow us to conclude that q ∈ C◦ . We can assume that pr∗ lies within a small gap arc(a, b), and that per∗ lies within a small gap arc(c, d). We will make use of these facts: (i) If a particular (axis-aligned) rectangle Di intersects Cq ’s boundary in a particular quadrant, then every smaller disc at q that intersects Di will also do so in that quadrant. (ii) If some Di intersects Cq ’s boundary in other than the first quadrant, then Di 6= Df . (iii) If the boundary of some Di intersects Cq at points u and v lying in the same quadrant, then arc(u, v) must lie within Di as well. 170 (iv) If the boundary of some Di intersects Cq at points u and v in opposing quadrants (i.e., quadrants 1 and 3, or 2 and 4), then Di must contain q, and the boundary of every smaller disc at q will intersect Di at points u0 and v 0 that have the same polar angle (with respect to q) as u and v respectively (Figure A.1). (v) If the boundary of some Di intersects Cq at points u and v that lie in adjacent quadrants (e.g., 1 and 2), then every smaller disc at q whose boundary intersects Di must to so at points u0 and v 0 , where u0 and v 0 lie within the (non-reflex) wedge formed by q, u, and v. u u0 q v0 Cq v Figure A.1: Fact (iv). Our proof will have a number of cases, based upon the quadrants containing points a, b, c, and d, and the rectangles of D associated with these points. For convenience, we will include a subscript with each of these points denoting its quadrant. For example, ‘b3 ’ refers to a point b lying within quadrant 3. Some of our arguments involve shrinking the original disc Cq . As the disc shrinks, its boundary will not contain any of the original gap endpoints (e.g., b3 ); but by fact (i), its boundary will contain some other point from the same rectangle, and in the same quadrant. For brevity, we will refer to this new boundary point by the same name and subscript (b3 ). There are 15 distinct cases (once the symmetric ones have been removed), based upon the quadrants containing the gap endpoints a, b, c, and d. Each of these cases has seven subcases, based upon which rectangles are associated with these endpoints (Table A.1). 171 Subcase 1 2 3 4 5 6 7 Properties δ(a) = δ(c) δ(a) = δ(d) δ(b) = δ(c) δ(b) = δ(d) δ(a) = δ(c) and δ(b) = δ(d) δ(a) = δ(d) and δ(b) = δ(c) a, b, c, and d lie in distinct rectangles of D Table A.1: Proof of Lemma 4.16, Subcases. Case a1 b1 c1 d3 Subcase 1 2,5,6 a1 b1 c2 d3 a1 b1 c2 d4 a1 b1 c3 d3 a1 b2 c2 d3 3,4 1 2,4,5,6 3 1,2 3,4 5,6 1 2,6 3,5 4 1,3 2,6 4 5 Proof steps One of (b1 , c1 , d3 ) or (d3 , a1 , b1 ) is a critical subset of Cq (to see this, consider rotating line pr per ccw around q until it reaches either b1 or d3 ). Shrink Cq until either (c1 , d3 ) is a critical subset; or until Df becomes outside-tangent to Cq , at which time (pr∗ , c1 , d3 ) is a −→ critical subset (d3 must remain to the right of − qp r ∗ , by fact (iv); c1 must remain to its left, since arc(c1 , d3 ) is not a large gap; and Df ∈ / {Dc1 , Dd3 }, by fact (ii)). (a1 , c1 , d3 ) is a critical subset. (b1 , c2 , d3 ) or (d3 , a1 , b1 ) (is a critical subset of Cq ). Shrink Cq until (pr∗ , c2 , d3 ) (is a critical subset of Cq ). (a1 , c2 , d3 ) (is a critical subset of Cq ). (b1 , c2 , d4 ). (a1 , c2 , d4 ). Shrink Cq until (c2 , d4 ) or (pr∗ , c2 , d4 ). (b1 , c3 , d3 ) or (d3 , a1 , b1 ). Shrink Cq until (c3 , d3 ) or (a1 , c3 ) or (pr∗ , c3 , d3 ). Shrink Cq until (c3 , d3 ) or (pr∗ , c3 , d3 ). (a1 , c3 , d3 ) or (c3 , a1 , b1 ). (d3 , a1 , b2 ). Shrink Cq until (c2 , d3 ) or (pr∗ , c2 , d3 ). (a1 , c2 , d3 ). Shrink Cq until (a1 , d3 ), (c2 , d3 ), (pr∗ , c2 , d3 ), or (pr∗ , d3 , a1 ). Table A.2: Proof of Lemma 4.16, Cases 1 . . . 5. In subcase 7, Cq must satisfy the support condition; hence we omit this subcase. The proofs for each of the 15 cases, and the remaining 6 subcases, are given in Tables A.2, A.3, and A.4. 172 Case a1 b2 c2 d4 a1 b2 c3 d1 a1 b2 c3 d3 a1 b2 c3 d4 a1 b3 c3 d1 a1 b3 c3 d3 a1 b3 c3 d4 a4 b2 c2 d3 Subcase 1 2,5,6 3,4 1,2 3,6 4 5 1,3 2,6 4 5 1 2,3,5,6 4 1 2 3 4 5 6 1 2,6 3 4,5 1,3 2,6 4 5 1,3 2 4,5 6 Proof steps (b2 , c2 , d4 ) or (d4 , a1 , b2 ). Shrink Cq until (c2 , d4 ) or (pr∗ , c2 , d4 ). (a1 , c2 , d4 ). (b2 , c3 , d1 ). Shrink Cq until (c3 , d1 ) or (pr∗ , c3 , d1 ) (c3 stays to the right of − −→ qp r ∗ by fact (v)). (a1 , c3 , d1 ) or (c3 , a1 , b2 ). Shrink Cq until (c3 , d1 ) or (pr∗ , c3 , d1 ). (d3 , a1 , b2 ). Shrink Cq until (pr∗ , b2 , d3 ). (a1 , c3 , d3 ) or (c3 , a1 , b2 ). Shrink Cq until (a1 , d3 ), (pr∗ , c3 , d3 ), or (pr∗ , d3 , a1 ). (b2 , c3 , d4 ) or (d4 , a1 , b2 ). Shrink Cq until (c3 , d4 ), (pr∗ , c3 , d4 ), or (pr∗ , b2 , c3 ). (a1 , c3 , d4 ) or (c3 , a1 , b2 ). (d1 , a1 , b3 ) or (b3 , c3 , d1 ). by fact (iii), (x ∈ arc(d1 , a1 ), b3 , c3 ). (x ∈ arc(b3 , c3 ), d1 , a1 ). (a1 , c3 , d1 ) or (c3 , a1 , b3 ). Shrink Cq until (pr∗ , c3 , d1 ). (x ∈ arc(d1 , a1 ), y ∈ arc(b3 , c3 )). (d3 , a1 , b3 ). Shrink Cq until (a1 , b3 ) or (pr∗ , b3 , d3 ). (a1 , x ∈ arc(b3 , c3 )) or a1 , c3 , d3 ). (a1 , x ∈ arc(b3 , d3 )). (d3 , a1 , b3 ). −→ qp Shrink Cq until (a1 , b3 ), (pr∗ , a1 , b3 ) (if b3 is right of − r ∗ ), or (pr∗ , b3 , d4 ) (otherwise). (c3 , a1 , b3 ) or (a1 , c3 , d4 ). Shrink Cq until (pr∗ , c3 , d4 ). (d3 , a4 , b2 ). Shrink Cq until (pr∗ , c2 , d3 ). −→ Shrink Cq until (pr∗ , c2 , d3 ) (if d3 is right of − qp r ∗ ) or (pr∗ , d3 , a4 ) (otherwise). Shrink Cq until (pr∗ , c2 , d3 ). Table A.3: Proof of Lemma 4.16, Cases 6 . . . 13. 173 Case a4 b2 c2 d4 a4 b2 c3 d3 Subcase 1 2 3 4 5 6 1,3 2 4 5 6 Proof steps (b2 , c2 , d4 ). Shrink Cq until (c2 , x ∈ arc(d4 , a4 )) or (pr∗ , c2 , d4 ). (x ∈ arc(b2 , c2 ), d4 , a4 ). (c2 , a4 , b2 ) or (a4 , c2 , d4 ). Shrink Cq until (pr∗ , c2 , d4 ). (x ∈ arc(d4 , a4 ), y ∈ arc(b2 , c2 )). (d3 , a4 , b2 ). −→ qp Shrink Cq until (pr∗ , c3 , d3 ) (if c3 is left of − r ∗ ) or (pr ∗ , b2 , c3 ) (otherwise). (c3 , a4 , b2 ). Shrink Cq until (pr∗ , b2 , c3 ), (pr∗ , c3 , d3 ), or (pr∗ , d3 , a4 ). Shrink Cq until (pr∗ , c3 , d3 ). Table A.4: Proof of Lemma 4.16, Cases 14, 15. 174
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Guarantees concerning geometric objects with imprecise...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Guarantees concerning geometric objects with imprecise points Sember, Jeffery 2011
pdf
Notice for Google Chrome users:
If you are having trouble viewing or searching the PDF with Google Chrome, please download it here instead.
If you are having trouble viewing or searching the PDF with Google Chrome, please download it here instead.
Page Metadata
Item Metadata
Title | Guarantees concerning geometric objects with imprecise points |
Creator |
Sember, Jeffery |
Publisher | University of British Columbia |
Date Issued | 2011 |
Description | Traditional geometric algorithms are often presented as if input imprecision does not exist, even though it is often unavoidable. We argue that in some cases, it may be desirable for geometric algorithms to treat this imprecision as an explicit component of the input, and to reflect this imprecision in the output. Starting with three problems from computational geometry whose inputs are planar point sets (Voronoi diagrams, convex hulls, and smallest bounding discs), we recast these as problems where each input point's location is imprecise, but known to lie within a particular region of uncertainty. Where algorithms to solve each of the original problems produce a single geometric object as output, the algorithms that we present typically produce either guaranteed or possible output objects. A guaranteed object represents qualities that can be guaranteed for every point set that is consistent with the uncertain regions, and a possible object represents qualities that exist for at least one such point set. By dealing with input imprecision explicitly, these guaranteed and possible objects can represent a more accurate view of what can be reliably inferred from the input data. |
Genre |
Thesis/Dissertation |
Type |
Text |
Language | eng |
Date Available | 2011-10-19 |
Provider | Vancouver : University of British Columbia Library |
Rights | Attribution-NonCommercial-NoDerivatives 4.0 International |
DOI | 10.14288/1.0052098 |
URI | http://hdl.handle.net/2429/38084 |
Degree |
Doctor of Philosophy - PhD |
Program |
Computer Science |
Affiliation |
Science, Faculty of Computer Science, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 2011-11 |
Campus |
UBCV |
Scholarly Level | Graduate |
Rights URI | http://creativecommons.org/licenses/by-nc-nd/4.0/ |
AggregatedSourceRepository | DSpace |
Download
- Media
- 24-ubc_2011_fall_sember_jeff.pdf [ 1.47MB ]
- Metadata
- JSON: 24-1.0052098.json
- JSON-LD: 24-1.0052098-ld.json
- RDF/XML (Pretty): 24-1.0052098-rdf.xml
- RDF/JSON: 24-1.0052098-rdf.json
- Turtle: 24-1.0052098-turtle.txt
- N-Triples: 24-1.0052098-rdf-ntriples.txt
- Original Record: 24-1.0052098-source.json
- Full Text
- 24-1.0052098-fulltext.txt
- Citation
- 24-1.0052098.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
data-media="{[{embed.selectedMedia}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.24.1-0052098/manifest