Implementation of Fortune’s Algorithms
Reports on a Course Programming Assignment
Computing the Voronoi Diagram[1]
Liu Zhenming(cs_lzm)[2]
1. Introduction
Voronoi Diagram has a large range of applications; typically through its assistance, we can answer the query of a nearest post office [1]. There is a commonly known algorithm, namely Fortune’s algorithm, which computes Voronoi Diagram in O(nlogn). In this project, I implemented this algorithm in C++ with some assistance of CGAL ([2]). Figure 1 shown below is the result of running the program on a randomly generated point set.
In this report, I will describe the hierarchy of the program, as well as the data structure and algorithms used and possible ways for further improvement. It is organized as follows: section 2 provides the flow of execution, section 3 describes the classes and functions involved in program, section 4 explores the implementation detail, and section 5-7 discussed the bugs, discoveries and possible improvements for this program, the appendix contains the source code of the whole program, which is also available in http://ihome.ust.hk/~cs_lzm/vord/source.htm. Throughout the whole report, I will use “I” and “we” alternatively, although the project is developed by individual.
Figure 1.1. Apply the program on a random generated point set
2. Program Overview
This program is implemented in an object-oriented approach with a least dependency of CGAL; meanwhile, I brought the concept of “event-driven” from windows programming ([3])., in particular, all the action performed during the execution are driven by events, which here are the point events and circle events. This technique makes the main flow of program clear. Program 2.1 illustrates a fragment in the main function. q, as an instance
247: while (crtEvent = q.ExtraMinY())
248: {
249: //cout << "crtEvent: " << endl;
250: //crtEvent->Output();
251: crtEvent->ActionPerform();
252: }
Program 2.1 A fragment from the main function.
The Event Driven Model makes the flow clear
of EventQueue, is the priority queue storing events; executing its member function returns a pointer to Event. Depending on the type of the event, Event::ActionPerform in line 251 will handle a site event or a circle event.
Besides the Event and EventQueue we see so far, a third major component throughout the program is the BeachLineHndl (sounds “Beach-Line Handler”). This stores the status of the beach-line, consisting of numerous arcs, in the Fortune’s algorithm ([1]). It also corresponds to the Status class in the implementation of segment intersection searching in [5]. Note that our major concern on the EventQueue and BeachLineHndl is the running time for each operation should be strictly bounded by O(logn). Section 4 will discuss how I implemented them in details.
A most effort is put in this program to ensure it is self-contained. Most basic components (e.g. lines, points) are defined and implemented. Besides we have to interact with Geomview using geometry objects in the CGAL, the only place we call CGAL functions is to compute the intersection of two lines, which could be replaced by codes available in any standard numerical analysis textbook. Consequently, the program is extremely easy to port to any object-oriented program
3. Class and their Functionality
Codes shown in Program 3.1 illustrate all classes and structures involved throughout the program. The entire definition of them is in line 81-224 (see appendix). The following part of this section will explain the features and functionality of these classes.
41: struct MyLine;
42: struct MyPoint;
43: struct ArcNode;
44:
45: class RangeComparator;
46: class BeachLineHndl;
47: class Event;
48: class EventQueue;
Program 3.1 Classes/Structs declared in the program
MyLine and MyPoint are the geometry classes in the planar graph. ArcNode is the structure representing the arcs in the beach line which stores the focus of the support parabola as well as the left and right adjacent sites. Meanwhile, all the arcs are stored in the linked list which represents the whole beach line in increasing x-order. Through this linked-list in the BeachLineHndl, the arc is able to decide its left boundary and right boundary in constant time. Readers may realize that a linked-list representation of beach line is not sufficient to support O(logn) searching operations. Basically we need one more balanced binary search tree associated with the linked list, in our implementation, we have chosen the map class in the STL([6], [7]) as the binary search tree. More discussion is available in section 4.
Note that this program does not maintain face relevant information for simplicity. This means the program only draws a Voronoi Diagram on the screen; no data-structure is maintained to store the diagram. Nevertheless, in this implementation, it is easy to find out the “owner” site for each segment, and hence not difficult add a face class in the future.
It is interesting to have a glance of the way we maintain the line class, which is called MyLine in our convention (See Program 3.2).
94: struct MyLine
95: {
96: MyLine(MyPoint _p1 = PointNULL, MyPoint _p2 = PointNULL): p1(_p1), p2(_p2)
97: {
98: terminate[0] = PointNULL;
99: terminate[1] = PointNULL;
100: refPt = PointNULL;
101: }
102:
103: void AddTerminate(MyPoint pt);
104: void Output() const;
105: MyPoint p1, p2;
106: MyPoint refPt;
107: MyPoint terminate[2];
108: };
Program 3.2 Definition of MyLine
The field p1 and p2 of type MyPoint indicate the two sites defining this segment (i.e. the segment is the perpendicular bisector of p1 and p2). Knowing that all the lines in the Voronoi Diagram are defined by two sites and introduced in the node-insertion, we only keep these two variables in the initialization. This avoids the massy mathematical or numerical stuff in the early stage. As the sweep line moves downward and the lines extend long enough and touch Voronoi vertices, we simply call AddTerminate to attach an end point for the corresponding lines. It provides an alternate, but easier way, to maintaining faces.
The RangeComparator stores either an arc in the beach line or a y-value. It serves as the key in the binary search tree to guide the search.
4. Data-Structure, Algorithms Used and the Running time analysis
4.1 Newton Method: In Fortune’s Algorithm, we need to decide the intersection of two parabolas precisely and frequently. This yields to solve quadratic equations. Having ever tried to apply standard formula, we realized unacceptable errors are introduced; therefore, we choose to use Newton method in solving quadratic equations. Since Newton method has a rapid convergence; normally returns a solution within 10 steps, we still consider it as a constant time operation.
4.2 Data Structure for BeachLineHndl: As discussed in section 2, we have maintained a linked list and a binary search tree. The linked list stores all the arcs in increasing x-order and is able to locate an arc’s neighbor quickly (which is impossible in the search-tree described in [1]); the binary search tree stores the node of the linked list and guides the search. This model is shown in Figure 4.1. In addition to the elegant property that a linked list owns, we can easily find out the desired data-structures in STL, getting rid of those headache, undocumented and unsafe code in the Internet.
4.3 Handling the Unclosed Point: After all the events are processed, there are segments(rays) that are unclosed (i.e. only one end point is attached). Our program is unable to determine the direction of these rays at this stage. We solve this problem by continue to move down the sweep line and trace the direction the rays move. This is finished in BeachLineHndl::BeachLineDone(Declared in line 193 and implemented in line 715. Till now we are able to tell their correct directions and draw them.
Figure 4.1 a Linked List and Balanced Search Tree that represents the beach line
5. Program Constraints (Bugs)
There are quite a few problems in this implementation. An obvious shortage is that the program only draws the Voronoi Diagram instead of computing the data-structure for it. Secondly the program will not work if any two points from the input have the same y-coordinate.
We encountered many numerical problems, and the solution for them is not elegant. For instance, sometimes, we realize that the left boundary even has a larger x-coordinate than the right boundary of a particular arc. We resolve this by simply using a large eps, and consider the close values be the same. Till the program is submitted, this method works fine; however, it could crash or output a wrong answer for the test data that the program does not like.
Lastly, knowing that not much resource can be released during the execution, I simply ignore to implement any destructor to collect the garbage. This works fine if this program is the only one to be executed, but surely will produce serious problems if it is involved as a part of other programs.
6. Possible Improvement
6.1 Less Quadratic Equations: Solving quadratic equations is expensive, especially using Newton’s method. It is possible to avoid some of them if we know we have already solved it before, and no change is made locally since then.
6.2 Separate File: Because of being unfamiliar with the compiling rule for cross-referencing for global variables, I have to use a large file to accomplish everything. It is highly undesirable since it produce difficulties in maintaining them.
6.3 Put More constrains in accessing variables: As the classes, variables are highly dependent on each other, changing status most like requires us to change other status (variables) elsewhere, this is difficult to avoid due to the complexity of the Fortune’s Algorithm. However, we can add more constrains on accessing the variables to avoid user operate on them by mistake
6.4 Numerical Errors: Till now I do not have any good ideas for it.
7. Conclusion
Loosely speaking, we have implemented the Fortune’s Algorithm successfully. However, the bugs stated above present it from being an industrial level product. Solving these problems requires the programmer a solid background in various aspects. This stimulates me pursuing further CS study. Lastly, it is really an unforgettable experience as I were working for it day after night.
Reference:
[1] M. deBerg, M. Van Kreveld, M. Overmars, O. Schwarzkopf, “Computational Geometry”, Springer 2000
[2] Computational Geometry Algorithms Library(CGAL) http://www.cgal.org
[3] Charles Petzold, “Programming Windows with C#”, Microsoft Press, 2002
[4] HKUST COMP 300B 2004 Spring Lecture Notes
http://www.cs.ust.hk/faculty/scheng/comp300b/notes/index.htm
[5] HKUST COMP 300B 2004 Spring Tutorial
http://www.cs.ust.hk/faculty/scheng/comp300b/tutorials/index.html
[6] SGI, Standard Template Library Programmer's Guide, http://www.sgi.com/tech/stl/
[7] / P.J. Plauger, “The C++ Standard Template Library”, Prentice Hall 2001
[8] Gerald Recktenwald, “Numerical Methods with MATLAB: Implementations and Applications”, Prentice Hall, 2000
Appendix the Source Code
[1] Report also available in http://ihome.ust.hk/~cs_lzm/vord
[2] Department of Computer Science, HKUST
