An Improved Algorithm for Computing the Volume of the Union of Cubes ∗
Pankaj K. Agarwal Department of Computer Science Duke University
Abstract Let C be a set of n axis-aligned cubes in R3 , and let U(C) denote the union of C. We present an algorithm that computes the volume of U(C) in time O(n polylog(n)). The previously best known algorithm takes O(n4/3 log2 n) time.
Categories and Subject Descriptors F.2.2 [Analysis of Algorithms and Problem Complexity]: Nonnumerical Algorithms and Problems—Computations on discrete structures, geometrical problems and computations
General Terms Algorithms, Theory.
Keywords Geometric data structures, Union of cubes.
1. INTRODUCTION Let C be a set of n axis-aligned cubes in R3 , and let U(C) denote the union of C. Computing the volume of U(C) efficiently is related to the well-known Klee’s measure problem. In 1977 Victor Klee [13] had presented an O(n log n) time algorithm for computing the union of n intervals in R1 and had asked whether his algorithm was optimal. An Ω(n log n) lower bound was proved by Fredmen and Weide [11]. Bentley [4] studied the two-dimensional version of Klee’s measure problem. When extended to computing the volume of the union of n d-dimensional axis-aligned boxes, the running time of Bentley’s algorithm is O(nd−1 log n). Later, van Leeuwen and Wood [12] improved the running time to O(n2 ) for d = 3. The problem lay dormant for a while until Overmars and Yap presented an algorithm with O(nd/2 log n) running time [15]. Recently, Chan [7] slightly improved the running time of their al∗ gorithm to nd/2 2O(log n) ; see also [8, 9] for some other related ∗ Supported
by NSF under grants CNS-05-40347, CCF-06 -35000, IIS07-13498, and CCF-09-40671, by ARO grants W911NF-07-1-0376 and W911NF-08-1-0452, by an NIH grant 1P50-GM-08183-01, by a DOE grant OEG-P200A070505, and by a grant from the U.S.–Israel Binational Science Foundation. Permission to make digital or hard copies of all or part of this work for personal or classroom use is granted without fee provided that copies are not made or distributed for profit or commercial advantage and that copies bear this notice and the full citation on the first page. To copy otherwise, to republish, to post on servers or to redistribute to lists, requires prior specific permission and/or a fee. SCG’10, June 13–16, 2010, Snowbird, Utah, USA. Copyright 2010 ACM 978-1-4503-0016-2/10/06 ...$5.00.
results. The lower bounds proved in [7, 16] suggest that improving the running time to O(nd/2−ε ), for any constant ε > 0, might be hard. See [1] for a brief history of the problem. Boissonnat et al. [5] proved that the combinatorial complexity of n axis-aligned cubes in Rd is Θ(n⌈d/2⌉ ), and it is Θ(n⌊d/2⌋ ) if all the cubes have the same size. Note that this bound is considerably better than the Θ(nd ) worst-case bound on the union of n axis-aligned boxes in Rd . This suggests that it might be easier to compute the volume of the union of n cubes. Indeed, the volume of the union of n axis-aligned unit cubes in R3 can be computed in O(n log n) time, by computing their union explicitly (which has linear complexity), but this will not lead to an efficient algorithm for cubes of different sizes. Agarwal et al. [2] presented an O(n4/3 log2 n) algorithm for computing the volume of the union of n axis-aligned cubes in R3 , by exploiting the special structure of cubes. Recently Bringmann [6] extended the result to higher dimensions and proposed an O(n(d+2)/3 ) algorithm for computing the volume of the union of n hypercubes in Rd . In this paper we significantly improve the result in [2]. T HEOREM 1.1. The volume of the union of a set of n axisaligned cubes in R3 can be computed in time O(n log4 n). Although there are some similarities between the new algorithm and [2], several new ideas and data structures are needed to improve the running time. The algorithm asserted in Theorem 1.1 is presented in the following four sections. Section 2 gives a highlevel description of the algorithm, and describes a data structure for maintaining the area of the union of a set of squares in R2 under certain assumptions on the sequence of updates. Section 3 describes the update procedure for the data structure, Section 4 presents auxiliary procedures needed by our algorithm in Section 3, and Section 5 presents secondary data structures needed in Section 4.
2.
THE GLOBAL STRUCTURE
We assume that the cubes of C are in general position. In particular, we assume that no plane support facets of two distinct cubes in C. Let z1 < · · · < z2n be the (distinct) z-coordinates of the vertices of cubes in C, sorted in increasing order. We sweep a horizontal plane Π in the (+z)-direction from −∞ to +∞, stopping at each zi . Let Π(t) denote the horizontal plane at z = t. For each 1 ≤ i < 2n, the cross-section U(C) ∩ Π(z) is the same for all z ∈ (zi , zi+1P ). Let αi denote the area of this cross-section. Then Vol U(C) = 2n−1 i=1 αi (zi+1 − zi ). We thus need to maintain αi as we sweep the horizontal plane. The intersection of Π(z) with U(C) is the union of a set S of squares that changes dynamically— a square is added to the intersection when Π sweeps through the bottom facet of its corresponding cube, and is removed from the intersection when Π sweeps through the top facet of its cube. Let
S1 and S2 be two squares that are inserted into S. If S1 is bigger than S2 and S1 is inserted after S2 , then S1 is deleted after S2 has been deleted (because S1 and S2 are intersections of cubes with the sweep line). We capture this observation by defining the notion of a size preserving update sequence:
A sequence of insertions and deletions is called size preserving if a square S is deleted only after all squares that were smaller than S and inserted before S have been deleted. We describe a dynamic data structure that maintains, in O(log4 n) amortized time (see Lemma 3.2), the area of the union of S, denoted by µ(S), for a size-preserving update sequence. This implies Theorem 1.1. 2v
µ(A ∪ B, ρ) = µ(B, ρ) + µ(A, B, ρ) .
(1)
Primary data structure. In our case, we know in advance the set S of all squares that will ever be inserted into S. S is the set of the xy-projections of cubes in C; |S| = n. Let 2 be the smallest axisparallel rectangle containing all squares in S. We build a binary space partition tree T on S, the so-called longest-sided kd-tree [10], that is similar to a kd-tree [3]. The structure of T will be static but the information stored at each node of T will change as the set S changes. + Each node v of T is associated with a rectangle 2v = [x− v , xv ]× [yv− , yv+ ] ⊆ 2 and a subset Sv ⊆ S of squares. For the root node u of T, 2u = 2 and Su = S. If the horizontal edge of 2v is longer + − − than its vertical edge (i.e., x+ v − xv ≥ yv − yv ), then we call v a -type node, otherwise we call v a ⊖-type node; see Figure 3. We classify a square S ∈ S that intersects 2v into one of the following mutually disjoint categories (see Figure 2): S is called a cover if 2v ⊆ S. S is an anchor if S contains exactly one edge of 2v — exactly one edge e of S intersects the interior of 2v . We set ϕ(S) := ϕ(e). If S contains a long edge of 2v , then S is called a long anchor, and a short anchor otherwise. S is a bar if 2v does not contain a vertex of S and two parallel edges of S—the ones parallel to a short edge of 2v —intersect the interior of 2v . If v is -type (resp. ⊖-type), then only the vertical (resp. horizontal) edges of S can intersect the horizontal (resp. vertical) edges of 2v . Finally, S is called a corner and a floater if 2v contains one vertex and at least two vertices of S, respectively, in its interior.
−
Figure 1. Long (dotted lines) and short (solid lines) squares in 2v , 2v (thick rectangle). Lightly shaded region is U(Lv ) ∩ 2v , and the hatched −
region is U(Gv ) ∩ 2v .
Here is a rough sketch of the overall data structure for maintaining µ(S). S is stored in a partition tree T, which is a variant of a kd-tree. Each node v of T is associated with a rectangle 2v . The squares of S intersecting 2v are partitioned into two sets Lv and Gv , the sets of long and short squares, respectively. A square is long if at most one of its edges intersects 2v , and short if at least
⊖
⊖
1111111111111 0000000000000 000 111 0000000000000 1111111111111 000 111 0000000000000 1111111111111 00 11 0000000000000 1111111111111 − 00 11 0000000000000 1111111111111 2v 0000000000000 1111111111111
gously. We will repeatedly use the following simple equality
bars long anchor
−
two edges intersects 2v . Let 2v = 2v \ U(Lv ), which is a rect−
corner
corner
angle; see Figure 1. The basic idea is to maintain Lv and 2v at v, to store Gv recursively at the descendents of v, and to compute
−
and even more challenging issue, is to maintain 2v and α(v) ˜ efficiently. A square can be long at too many nodes of T, so a long square is stored only at the highest nodes at which it is long (anal−
ogous to how intervals are stored in a segment tree). But then 2v cannot be not maintained explicitly at v, making the maintenance of α(v) ˜ hard. We circumvent this problem by maintaining an ap−
proximation ⊞v of 2v , maintaining α(v) = Area U(Gv ) ∩ ⊞v , pushing certain information related to Lv to the descendents of v in a lazy manner, and maintaining secondary structures to compute α(v) ˜ from α(v) quickly whenever needed. A clever charging scheme is used to bound the time spent in maintaining this information. The charging scheme exploits, in a crucial though subtle way, the special (obvious) property that the life-time of a square of side length h (regarding the z-direction as “time”) is also h time units. For a horizontal (resp. vertical) edge e lying on the line y = a (resp. x = a), we refer to a as the value of e and denote it by ϕ(e). For two sets A and B of squares and a rectangle ρ, let U(A, B, ρ) denote (U(A) ∩ ρ) \ U(B), and µ(A, B, ρ) denote Area U(A, B, ρ). For brevity, we set U(A, ρ) := U(A, ∅, ρ) and U(A, B) := U(A, B, R2 ). We define µ(A, ρ) and µ(A, B) analo-
short anchor
floaters
Figure 2. Partition of S into various categories.
We refer to bars, corners, and floaters as short squares (at least two edges of a short square intersect 2v ), and to anchors and covers as long squares. We call a square S active at v if it is not a cover or a long anchor at any proper ancestor of v. All short squares at v are active, and an active anchor at v is allowed to be an active short anchor at p(v), the parent of v. Let Sv ⊆ S be the subset of active squares at v. For a square S ∈ S, let Ξ(S) be the set of lowest nodes of T at which S is active, and let N(S) be the set of (proper) ancestors of nodes in Ξ(S). S is active at all nodes in N(S) ∪ Ξ(S), and it is stored at these nodes. Similarly, we call an orthogonal segment e intersecting 2v short at v if 2v contains an endpoint of e, and long otherwise. We call e active at v if there is no proper ancestor u of v at which e is long and parallel to a long edge of 2w ; this definition is similar to a square being active at a node v. Finally, we define Ξ(e) and N(e) for an edge e as for a square. If there are no short squares at v, then v is a leaf. Otherwise, v is an interior node, and we split 2v into two rectangles by splitting its longer edges, as follows. Suppose v is -type. Let Xv be the set of x-coordinates of the vertices of short squares at v that lie in the + interval [x− v , xv ], i.e., the set of x-coordinates of the vertices of: (i) S that lie inside 2v , and (ii) of bars at v. Let χv be the median of Xv . We split 2v into two rectangles by drawing the vertical ⊖
−
α(v) ˜ = Area U(Gv ) ∩ 2v from the information stored at the children of v. There are, however, several challenges in accomplishing √ this: First, a square may appear short at Ω( n) nodes in a kd-tree, so how do we ensure that a square appears as a short square only at polylog(n) nodes of T? Exploiting that S is a set of squares, we describe a variant of kd-tree that guarantees this property. Second,
line ℓv : x = χv . The resulting left (resp. right) subrectangle is associated with the left (resp. right) child of v. If v is ⊖-type, then we split 2v by the horizontal line ℓv : y = χv , where χv is the median of the y-coordinates of the vertices of short squares at v that lie in the interval [yv− , yv+ ]; see Figure 3. We associate the bottom (resp. top) subrectangle with the left (resp. right) child of v. We call ℓv the separator line at v. We define a subtree T(v) of T rooted at v, as follows: we recursively visit the subtree rooted at v. Suppose we are at a node w. If w is a leaf of T or a node of different type than v, we stop. Otherwise, we recursively visit the children of v. T(v) is the set of nodes visited by this procedure. Let D(v) denote the set of leaves of T(v). All interior nodes of T(v) are of the same type as v. See Figure 3.
path in ΠB , or S appears as a corner at w. Hence, the number of paths in ΠA is O(∆), and S appears as an active short anchor at O(∆2 ) nodes. If S appears as an active cover at a node v, then S appears as an active short anchor at p(v), implying that S appears as an active cover at O(∆2 ) nodes. Finally, if S appears as a long anchor at v, then either p(v) is the leaf of a path in ΠA ∪ ΠB or S appears as a corner at p(v). Since there are O(∆) such nodes, S appears an active long anchor at O(∆) nodes. Next, we prove that ∆ = O(log n), which will complete the proof of (ii). c r
d
C Da B b c A
H h
F g
f G
E
K i e
jJ I
f
e
g
d
2r
g b
a
a
a c f
d
b
b c
A B C
i
d
f D
E
g
h
j
K
FG HI J
Figure 3. A recursive partition of 2 and the corresponding T; lower-case letters denote the separator lines, and upper-case letters denote the cells associated with the leaves of T; shaded region denotes T(rt).
Each square is a floater at the root of T. Let S ∈ S be a square, and let v be a node such that S is active at v and the separating line ℓv intersects the interior of S. If S ⊂ int 2v , S splits into two floaters. Assume that S intersects ∂2v . If S is a floater at v, then either S splits into two corner squares, or into a bar and a floater. If S is a bar at v, it splits into two (active) anchors. (Note that since ℓv splits the long edges of 2v , a bar is never split into two bars.) If S is a corner at v, then it splits into an (active) anchor and a corner. Finally, if S is an active short anchor at v, it splits into a cover and an anchor. See Figure 4. The following lemma will be used to bound the size of the data structure. L EMMA 2.1. (i) T has O(n log n) leaves, and the depth of T is O(log n). (ii) A square appears as a floater, a corner, or an active long anchor at O(log n) nodes, and as a bar, an active short anchor, or an active cover at O(log2 n) nodes. Proof: Let ∆ be the depth of T, and let S be a square in S. If S is a floater or a corner at a node v, then 2v contains at least one vertex of S. Since there are at most four nodes v of T at any given level for which 2v contains a vertex of S, S appears as a floater or a corner in at most 4∆ nodes. Next, if S appears as a bar at a node v, then S appears at p(v) as a bar or a floater and does not appear as a bar at the sibling of v. Therefore the nodes at which S appears as a bar can be partitioned into a family ΠB (S) of paths, each of length at most ∆, such that S appears as a floater at the parent of the root of each path. See Figure 4. Since S appears as a floater in at most 4∆ nodes, S appears as a bar in at most 4∆2 nodes. If S appears as an active short anchor at v, then there are two cases: (i) S appears as a bar at p(v), S intersects ℓp(v) , and S appears as an anchor at the sibling of v, or (ii) S appears as a corner or an active short anchor at p(v) and does not appear as a short anchor at the sibling of v. We can thus partition the nodes at which S appears as an active short anchor into a family ΠA of paths each of length at most ∆ so that either the parent w of the root of each path is a leaf of a
Figure 4. Trace of a square S through T. Left: S appears as a corner or floater along thick paths; paths in ΠB (S) (resp. ΠA (S)) are denoted by dashed (resp. solid) lines; S appears as a cover (resp. long anchor) at light (resp. dark) squares, and as a short anchor (resp. bar) at hollow (resp. filled) circles. Right: Separator lines inside 2r corresponding to the marked nodes on the left.
Let L(m, k) denote the maximum number of leaves in a subtree rooted at a node v so that Sv contains m bars and 2v contains k vertices of S. Let ∆(m, k) denote the maximum depth of such a subtree. If m + k ≤ 1, then ∆(m, k) ≤ 1 and L(m, k) ≤ 2. For m + k ≥ 1, let w1 , w2 be the two children of v. For i = 1, 2, let mi , ki denote the number of bars at wi and the number of vertices that lie in 2wi . Then L(m, k)
≤ L(m1 , k1 ) + L(m2 , k2 ),
∆(m, k)
≤ 1 + max{∆(m1 , k1 ), ∆(m2 , k2 )}.
A bar at v contributes at most one bar at a child of v, and a floater at v contributes at most one bar at a child of v. At least two of the vertices of a floater at v lie in 2v , so the number of floaters at v is at most k/2, thereby implying that m1 + m2 ≤ m + k/2 and k1 +k2 ≤ k. Since the separator line chooses the median of 2m+k values, namely the (x- or y-) coordinates of the vertices of bars and of the vertices that lie in 2v , 2mi + ki ≤ (2m + k)/2 for i = 1, 2. The solution to the above recurrence is ∆(m, k)
≤
A1 log2 (2m + k),
L(m, k)
≤
A2 (m + (k/2) log2 (2m + k)),
where A1 , A2 ≥ 1 are constants. Since there are no bars at the root of T and S has at most 4n vertices, the depth of T is O(log n) and it has O(n log n) leaves, as claimed. 2 Information stored at a node v. Let S ⊆ S be the current set of squares. Let Sv = Sv ∩ S be the subset of current squares that are active at v; T stores Sv at v. Let Bv (resp. Cv , Fv ) be the subset of squares in Sv that are bar (resp. corner, floater) squares at v. Set Gv := Fv ∪ Bv ∪ Cv ⊆ Sv to be the set of short squares of Sv . Let Xv (resp. LAv , SAv ) be the subset of active covers (resp. long anchors, short anchors) at v. Set Lv := LAS v ∪ S A v ∪ Xv ⊆ Sv to be the set of long squares of Sv , and L∗v = u Lv where the union is taken over all ancestors of v including v itself. The following lemma is straightforward:
=
µ(L∗v , 2v ) + µ(Gv , L∗v , 2v ) −
−
Area(2v ) − Area(2v ) + µ(Gv , 2v ) . −
See Figure 5 (left). Ideally, we would like to maintain 2v and −
µ(Gv , 2v ) at each node v, from which µ(S, 2v ) can be computed. It is, however, too expensive to maintain them explicitly at each −
⊖
⊖
node because the insertion or deletion of a square may change 2 at several nodes. We therefore maintain another window ⊞v called − exposed window, an approximation of 2v defined below, and maintain µ(Gv , ⊞v ). + If v is -type (resp. ⊖-type), then we define L B− v , L B v to be the + bottom and top (resp. left and right) edges of 2v , and S B− v , S B v to be the left and right (resp. bottom and top) edges of 2v —the former two are the long edges of 2v and the latter two are the short edges + − + − of 2v . Let lb− v , lbv , sbv , and sbv denote the values of L B v (i.e., − + − + lb− = ϕ(L B )), L B , S B , and S B , respectively. 2 v is thus v v v v v − + − + + − + [sb− v , sbv ] × [lbv , lbv ] if v is -type, and [lbv , lbv ] × [sbv , sbv ] otherwise. We divide active long and short anchors at v into two + sets each: LA− v (resp. L A v ) is the subset of L A v that contains the − + + edge L Bv (resp. L Bv ), and SA− v (resp. S A v ) is the subset of S A v + that contains the edge S B− (resp. S B ). To maintain ⊞v and to v v −
compute 2v , we maintain the following information at v: I. The sets Fv , Cv , Bv , LAv , SAv , and the size |Xv |. − II. The sequence La− v = hϕ(S) | S ∈ L A v i in decreasing + order, and the sequence Lav = hϕ(S) | S ∈ LA+ v i in in+ creasing order. For technical reasons, we add lb− v (resp. lbv ) − + at the end of Lav (resp. Lav ). − III. The sequence Sa− v = hϕ(S) | S ∈ S A v i in decreasing order, and the sequence Sa+ = hϕ(S) | S ∈ SA+ v v i in increas+ ing order. For technical reasons, we add sb− v (resp. sbv ) at − + the end of Sav (resp. Sav ). IV. Fv : the edges of squares in Fv parallel to the long edges of 2v sorted in increasing order of their values. −
+ Anchored window. 2v is also defined by four edges L A− v , L Av , −
+ − + − + S A− v , S A v , and let lav , lav , sav , sav be their values. If 2 = ∅, then these values are not well defined. We define them directly −
−
so that they are well defined even if 2 = ∅, and view 2v as the intersection of four halfplanes defined by the lines supporting ˜ these edges. We describe la+ ϕ, and let v : Set λv = minϕ∈La+ v + ˆ v = la ˆ v = sa+ λ if p(v) and v are of the same type and λ p(v)
p(v)
˜ ˆ otherwise. Then la+ v is defined to be the minimum of λv and λv but + − + its value cannot be less than lb− (to ensure that la ∈ [lb , v v v lbv ]), − ˆ ˜ i.e., la+ = max{lb , min{ λ , λ }}. Values of the other edges v v v v −
−
of 2v are described similarly. These values, and thus 2v , can be computed by traversing the path in T from the root to v. Exposed window. The exposed window at v, ⊞v ⊆ 2v , is also + − + composed of four edges, called virtual walls, L E− v , L Ev , S Ev , S Ev − + − + and whose values are lev , lev , sev , sev , respectively. See Figure 5. For each node, the algorithm maintains these edges and the quantity α(v) := µ(Gv , ⊞v ) = µ(Fv ∪ Bv ∪ Cv , ⊞v ),
(2)
−
If ⊞v is tight, then α(v) = µ(Gv , 2v ), as desired. As mentioned above, it is too expensive to ensure that ⊞v is tight at all nodes of T. The algorithm lets ⊞v be loose at many nodes, and updates them carefully in a lazy manner so that the following two invariants are maintained at every node v ∈ T; see Figure 5 (ii): (I1) If v is -type (resp. ⊖-type), then the y -coordinate (resp. x-coordinate) of none of the vertices of a square in Fv lies − + + between le− v and lav or between lev and lav . (I2) If v is -type (resp. ⊖-type), then the x-coordinate (resp. y coordinate) of none of the vertices of a square in Fv ∪ Bv lies − + + between se− v and sav or between sev and sav . Figure 5 (iii) illustrates an example of an exposed window that violates (I1) and (I2). Intuitively, (I1) and (I2) ensure that no edge + + of Fv ∪ Bv that is parallel to L E+ v lies between L E v and L A v , though an edge of Cv may lie there; the same holds for the other three virtual walls of ⊞v . These invariants enable us to compute −
µ(Fv ∪ Bv , ρ), for any ρ ⊆ ⊞v ⊕ 2v , by reducing it to a 1D problem. The secondary data structures described below will facilitate the computation of µ(Gv , ρ) (provided I1 and I2 are maintained), which in turn is used to update α(v) as ⊞v is made tight. Secondary structures at v. For a square S, with a slight abuse of notation, we use Sv− (resp. Sv+ ) to denote the interval in R1 formed by the intersection of S ∩ 2v and the line containing the edge L E− v − − (resp. L E+ v ) of ⊞ v . (See Figure 6.) Let Iv = {Sv | S ∈ Bv ∪ Fv }, + − + and let I+ v = {Sv | S ∈ Bv ∪ Fv }. For an interval I ∈ Iv ∪ Iv , − + q − + let I be the rectangle I × [yv , yv ] if v is -type and [xv , xv ] × I q − + q + if v is ⊖-type. Set P− v = {I | I ∈ Iv } and Pv = {I | I ∈ Iv }. The following two secondary data structures are maintained at v. Wall structure. Store I− v into a dynamic segment tree that supports the following two operations in O(log n) time [3]: (i) insertion or deletion of an interval to I− v , and (ii) for a query interval δ, + return U L− (v, δ) = |U(I− v )∩δ|. Similarly, store Iv into a segment + + tree so that I+ can be updated and U L (v, δ) = |U(I v v ) ∩ δ|, for a query interval δ, can be computed in O(log n) time. Corner structure. We process Cv and I− v into a data structure that supports the following operations: (i) insertion or deletion of a square to Cv and an interval to I− v ; (ii) for a query rectangle ρ ⊆ 2v , compute C A− (v, ρ) = µ(Cv , P− v , ρ) . Similarly, main+ tain Cv and I+ v into a data structure to compute C A (v, ρ) = µ(Cv , P+ , ρ) for a query rectangle ρ. The corner data structure v is described in Section 5 and used in Section 4. By Lemma 5.1, a query is answered in O(log n) time and an update can be performed in O(log2 n) time. ⊖
µ(S, 2v ) =
−
For a leaf v and for a node v with ⊞v = ∅, α(v) = 0. If ⊞v = 2v , then ⊞v is called tight (see Figure 5 (i)); otherwise it is called loose.
⊖
−
We call 2 := 2v \ U(Lv ∗ ) an anchored window, which is a (possibly empty) rectangle lying inside 2v . By (1),
i.e., the area of the union of short squares (at v) lying inside ⊞v . If + − + se− v ≥ sev or lev ≥ lev , then ⊞ v = ∅. For the root rt of T, the algorithm ensures that ⊞rt = 2rt , implying that α(rt) = µ(S).
⊖
L EMMA 2.2. Let w be a descendent of v. Then Fw ⊆ Fv and Gw ⊆ Gv . Moreover, if w is a node of T(v) of the same type as v, then Bw ⊆ Bv , otherwise Bw ⊆ Fv .
3.
UPDATE PROCEDURES
This section describes the procedures for updating T and the secondary structures stored at its various nodes, as a square is inserted or deleted. The update procedures maintain (I1) and (I2) and use the following auxiliary procedures. + − + A DJUST WALL(v, W L, ξ). Shift W L ∈ {L E− v , L Ev , S Ev , S Ev } of the exposed window ⊞v at v from its current position to ξ. The value of ξ is such that the invariants (I1) and (I2) are not violated by the shift. The procedure returns the updated value of α(v).
L E+ v 00000000000 11111111111 S E− v
00000000000 11111111111 11111111111 00000000000 00000000000 11111111111 00000000000 11111111111 00000000000 11111111111 00000000000 11111111111 00000000000 11111111111 L E− v 00000000000 11111111111 00000000000 11111111111
+ (se+ v , lev )
se− v 00000000000 11111111111 11111111111 00000000000 00000000000 11111111111 00000000000 11111111111 00000000000 11111111111 00000000000 11111111111 00000000000 11111111111
S E+ v
le− v la− v
− (se− v , lev )
sa− v
Figure 5. Exposed window ⊞v (lightly shaded rectangles) and its virtual walls, hatched region is U(Gv , ⊞v ). Left: ⊞v is tight; middle: ⊞v satisfies invariants (I1) and (I2); right: ⊞v violates both (I1) and (I2) — shaded strips contain vertices of floaters.
These procedures are described in Section 4. A DJUST WALL and T IGHTEN W INDOW procedures take O(log2 n) time, and P USH U P takes O(1) time. Assuming that we have these procedures at our disposal, we describe the procedure for inserting a square; the deletion procedure is symmetric. Inserting a square. Let S be a square that we wish to insert into S. We update the information stored at various nodes of T, including α(·). The new value of α(rt) is precisely the new value of µ(S). We first visit the nodes of N(S) and their children in a top-down manner, starting from the root, and call the function T IGHTEN W INDOW at each node. This ensures that the exposed window is tight at all these nodes before S is actually inserted. Next, for each node v ∈ N(S) ∪ Ξ(S), we insert S at v, as described below, in a top-down manner. Finally, we visit the nodes v ∈ N(S) ∪ Ξ(S) in a bottom-up manner to update the value of α(v). For v ∈ Ξ(S), α(v) is updated directly as described below, and for v ∈ N(S), it is updated using the procedure P USH U P(v); the exposed windows at v and its children will be tight before P USH U P(v) is invoked at v; it ensures that α(v) is correctly computed. The processing of S at a node v ∈ N(S) ∪ Ξ(S) depends on the status of S at v. Recall that if S is a cover or a long anchor at v, then v ∈ Ξ(S); if S is a floater, corner, or bar, then v ∈ N(S); if S is an active short anchor and v is a leaf (resp. an interior node) of T, then v ∈ Ξ(S) (resp. v ∈ N(S)). If S is a short square or a short anchor at v, then we simply insert S at v and update the information (I)-(IV) and the secondary structures at v in a straightforward manner in O(log2 n) time. Inserting an active long anchor or cover at v is more involved, as it may violate (I1) and (I2) at many of its descendents. Additional processing is needed at these descendents to restore (I1) and (I2). We describe these two cases in detail. la+ u le+ u
le+ v
le+ v ⊞v
⊞u
⊞w
C B
⊞u
S la+ u
⊖
some of the descendents of v to ensure that (I1) and (I2) are not violated after the insertion of S. We refer to this procedure as P USH D OWN; it also updates the value of α(·). For the sake of concreteness, we assume that S ∈ LA+ v and that v is -type; this implies −
+ that L A+ v , L E v are the top edges of 2v , ⊞ v . T IGHTEN W INDOW was called at v in the beginning of the insertion procedure, so ⊞v was tight before the insertion of S. If ϕ(S) ≥ la+ v , then ⊞ v remains tight even after the insertion of S and we do nothing. If + + ϕ(S) < la+ v then L A v lowers to ϕ(S). We lower L E v of ⊞ v to ϕ(S), by sweeping from its current position to ϕ(S) and stopping at all values in Fv , which we refer to as events. An event due to an edge e ∈ Fv indicates that (I1) or (I2) is violated at v (and some of its descendents) because of e. For technical reasons, we lower + L A+ v , L E v simultaneously through the sweep process. Let ξ denote the value of le+ v at the current event and η at the last event. Initially, we set ξ to the original value of le+ v (before the insertion of S) and + η = +∞. If L E+ v reaches ϕ(S), set α(v):=A DJUST WALL (v, L E v , + + + ϕ(S)) and lev = lav = ϕ(S), and stop. If L Ev reaches an edge e ∈ Fv , where e is an edge of a floater σ ∈ Fv , set η to ξ and ξ to ϕ(e). Let Z(v, e) be the subtree consisting of the descendents w of v that lie in N(e), at which σ is floater or bar, and for which its top virtual wall if above e. Z(v, e) can be computed in O(log2 n) time by traversing T in a top-down manner, starting from v. Let w be a node in Z(v, e). The insertion of S causes the top edge −
of 2w to pass below γ and thus the top virtual wall of ⊞w , denoted by W L, needs to be aligned with e to maintain (I1) and + (I2); W L is S E+ w (resp. L E v ) if v is ⊖-type (resp. -type). We set α(w):=A DJUST WALL(w, W L, ξ). If w is -type, we also update the secondary data structures at w: we insert the interval σv+ into + I+ w if e is the top edge of σ (L E w starts intersecting σ as we lower + it beyond e) and delete it from I+ w if e is the bottom edge (L E w no longer intersects σ as we lower it beyond e), and we update the wall and corner data structures at w. See Figure 6. For each node w ∈ Z(v, e), we say that S encountered e at node w. The procedure spends O(log2 n) time in computing Z(v, e) and O(log2 n) time at each node w ∈ Z(v, e); we charge the former to the “encounter” of e at v. For S ∈ LA− v and for a ⊖-type node, the P USHP D OWN procedure is symmetric. The total running time is O((1 + w κw ) log2 n), where κw is the number of edges encountered by S at w. ⊖
−
becomes the same as 2v . The procedure returns the updated value of α(v). P USH U P(v). Let w, z be the children of v. Assuming ⊞v , ⊞w , ⊞z are tight, the procedure computes the value of α(v) from α(w) and α(z) and returns it.
⊖
T IGHTEN W INDOW(v). Tighten the exposed window ⊞v so that it
a A
⊞w
Figure 6. Insertion of a long anchor. Solid (dashed) intervals belong to + I+ u (Iv ). Left: before the insertion of S. Right: after the insertion of S; + P USH D OWN procedure lowers L E+ u to the bottom edge of S, L E v to the top edge of B, and L E+ to the top edge of C. w + Long anchor. If S ∈ LA+ v , we insert ϕ(S) into Lav , and into − Lav otherwise. Next, we update the exposed windows at v and
x− v
b
c B Iv
C
S
u a
d D
v c
b Z d
S E+ u A B C D
Figure 7. Insertion of a cover square and the subtree Z. Left: shaded region is 2v , dotted line is S E+ v.
Cover. We increment the value of |Xv |. If |Xv | = 1, we update ⊞v and further processing is needed at v and its descendents to
⊖
⊖
ensure (I1) and (I2) and to update the secondary structures. Let u = p(v). Suppose u is -type and v is the right child of u, in which case S ∈ SA+ u ; see Figure 7. We move the right virtual wall + of ⊞v (which is L E+ v if v is ⊖-type and S E v otherwise) from its − current position, say, ξ, to xv , the value of the left boundary of 2v . Let Iv = [x− v , ξ], and let Z ⊆ T(u) be the set of nodes z such + that [x− z , xz ] ∩ Iv 6= ∅ and Fz ∪ Bz 6= ∅. Z forms a subtree of T(u) rooted at v; if v is ⊖-type, then Z = {v}. The total time spent in computing Z is O((1 + βv ) log n) where βv is the number of edges in Fv ∪ Bv whose values lie in the interval Iv . To ensure (I1) and (I2), we call the P USH D OWN procedure at the leaves w of Z that are of -type to move L E+ w to the left edge of 2w . At other nodes w of Z, we set α(w) = 0 and S E+ w to the left edge of 2w . The insertion procedure for other cases—when v is the left child of u or when u is ⊖-type—is P similar. The total time spent in inserting a cover at v is O((1 + w κw ) log2 n + βv log n), where w is a descendent of a leaf of Z and κw is the number of edges encountered by the P USH D OWN procedure at w. A square can be deleted using a similar procedure. The correctness of the update procedures follows from the following lemma. The proof of the lemma, tedious but not difficult, is omitted from this extended abstract because of lack of space.
+ implying that la+ v , and thus law , remain less than ϕ(γ) until σ is deleted. Since σ ∈ Fw and invariants (I1) and (I2) are maintained, le+ w ≤ ϕ(γ) until γ is deleted. Hence, γ will not be encountered again at w while lowering L E+ w . Similarly, γ will be encountered at w at most once by the insertion procedure while raising L E− w. Similar arguments show that γ is charged at w at most twice by the deletion procedure or when γ is parallel to the short edges of 2w . 2
Remark. The proof of the above lemma crucially uses the observations that if the P USH D OWN procedure at v because of the insertion of a long anchor S sweeps through the edge γ of a floater square σ then S is bigger than σ, and that the size-preserving order ensures that an edge is encountered only O(1) times for a fixed node. These claims are not true if P USH D OWN is called to move a short virtual edge. That is why short and long anchors (and edges of 2v ) are treated differently and auxiliary data structures are stored only with respect to long edges.
4.
AUXILIARY PROCEDURES
If α(v) 6= 0 and ⊞v , ⊞w , ⊞z are tight, where w and z are the children of v, then
L EMMA 3.1. Invariants (I1)–(I2) are maintained at each node of T after inserting a square into S or deleting from it.
α(v) = α(w) + α(z) + Area(⊞v \ (⊞w ∪ ⊞z )).
Since (I1) and (I2) hold at all times, the auxiliary procedures compute the value of α(·) correctly, which in turn implies that the algorithm correctly maintains α(v) at each node v of T at all times. In particular α(rt) = µ(S) at all times. The next lemma bounds the amortized update time.
Using this observation, P USH U P(v) can be implemented in O(1) time. T IGHTEN W INDOW can be implemented by invoking A D JUST WALL at most four times. We thus describe A DJUST WALL , which uses the secondary structures (wall and corner structures) stored at the nodes of T and relies on invariants (I1) and (I2).
L EMMA 3.2. We can maintain S into a data structure so that µ(S) can be updated in O(log4 n) amortized time, as we insert or delete a square, provided that the sequence of updates is size preserving. Proof: The argument in the update procedure implies inP that the 4 2 sertion/deletion time of a square S is O(log n + κ log n + v v P w βw log n), where κv is the number of edges encountered at v by S during the update procedure, w is a node at which S is a cover, and βw is as defined above. We charge O(log4 n) time to the insertion of S. We prove that the amortized cost of the second term is also O(log4 n); a similar argument bounds the amortized cost of the third term by O(log4 n). We charge log2 n units to each edge encountered by S during its insertion or deletion. We claim that an edge γ of a square σ ∈ S is charged at most twice at a node w of T by the insertion and deletion procedures each. Since w ∈ N(γ) and σ ∈ Fv ∪ Bv , γ is charged a total of O(log2 n) times. Hence γ is charged a total of O(log4 n) units, which can be charged to the insertion of σ. This implies that the amortized update time of a square is O(log4 n). To prove the above claim, suppose γ is parallel to the long edges of 2w ; the other case is symmetric. First consider the insertion procedure. Suppose γ was encountered by S at w while lowering le+ w (actually, γ ∈ Fw in this case). Recall that γ is encountered at w because the P USH D OWN procedure swept through γ at an + ancestor v of w while lowering L E+ v , L A v , where γ ∈ Fv and S ∈ LAv ∪ Xv . Since σ ∈ Fv , one of the edges of σ lies inside 2v while S contains a long edge of 2v , which implies that |X(σ)| < − (le+ v − lev ) ≤ |X(S)|, i.e., S is bigger than σ. (Note that this claim might not be true if S is a short anchor at 2v and thus only contains a short edge of 2v .) Moreover, σ is inserted into S before S. We have assumed the update sequence to be size preserving (see Section 2), therefore σ will be deleted before S is deleted,
000000 111111 111111 000000 000000 111111 000000 111111
e
f
L E+ v
ξ ρ
11111 00000 1111 0000 00000 11111 0000 00 00000LE+v 1111 11111 000011 1111 ⊞v
⊞v 2v
2v
+ Figure 8. A DJUST WALL(v, L E+ v ,ξ): current ⊞v and after adjusting L E v (lightly shaded rectangles). Left: Segment e and rectangle ρ are compatible + with L E+ v ; segment f is not compatible with L E v . Right: area of the hatched region is the change in α(v).
A DJUST WALL(v, W L, ξ). This procedure shifts the virtual wall W L of ⊞v from its current position to ξ, assuming that the value of ξ is such that the invariants (I1) and (I2) are not violated by the shift. The procedure returns the updated value of α(v). See Figure 8. Let R ⊆ 2v be the rectangle obtained by sweeping W L from its current position to ξ, then α(v) changes by µ(Gv , R). We call an orthogonal segment e ⊆ 2w compatible with a virtual wall W L of ⊞v parallel to e if the strip We bounded by the lines supporting e and W L does not contain a vertex of any square in Fv ∪ Bv . That is, no square of Fv ∪ Bv is contained in We and the two boundary lines of We intersect the same set of squares in Fv ∪ Bv . We call a rectangle ρ ⊆ 2v compatible with W L if the two edges of ρ that are parallel to W L are compatible with e. See Figure 8. For a square ρ ⊆ 2v that is compatible with one of the virtual walls of ⊞v , we describe a procedure that computes µ(Gv , ρ). By (I1) and (I2), R is compatible with W L, so µ(Gv , R) can be computed using this procedure. First, consider the case when ρ is compatible with L E+ v . In this case, ρ intersects a square S ∈ Fv ∪ Bv q if the interval I := Sv+ lies in I+ v . Furthermore, S ∩ ρ = I ∩ ρ, S
5.
be regarded as a bar. Let ρ = xρ × yρ . By (1), µ(Gv , ρ)
= µ(Fv ∪ Bv , ρ) + µ(Cv , Fv ∪ Bv , ρ)
= =
+ µ(P+ v , ρ) + µ(Cv , Pv , ρ) U L+ (v, xρ )|yρ | + C A+ (v, ρ).
Similarly if ρ is compatible with L E− v , then µ(Gv , ρ) = U L− (v, xρ )|yρ | + C A− (v, ρ). Using the wall and corner structures stored at v, µ(Gv , ρ) can be computed in O(log n) time. ρ2
v0
ρ
v4
111111111111 00000000 0000 00000000 11111111 0000 1111 00000000 11111111 0000 000000001111 v111111111 v2 v3
⊞v
I
se+ v
Figure 9. Computing µ(Gv , ρ) for a rectangle ρ compatible with S E+ v; thick rectangle is ρ and dashed lines are separator lines.
⊖
Next, suppose ρ is compatible with S E+ v and v is -type. Since we do not store secondary structure for short edges of 2v , answering the query is more involved in this case. By (I2), no vertical edge of a square in Fv ∪ Bv intersects ρ. Hence, if a square in Bv intersects ρ, it contains ρ and µ(Gv , ρ) = Area(ρ). We can check this condition in O(log n) time using the wall data structure. So assume that ρ does not intersect any square in Bv . Let v0 (resp. vk ) be the node in D(v) (cf. Section 2) such that 2v0 (resp. 2vk ) contains the left (resp. right) edge of ρ. Let I be the x-projection of ρ \ (2v0 ∪ 2vk ), and let v1 , . . . , vk−1 be the nodes in Ξ(I); each vi lies S in T(v). See Figure 9. For 1 ≤ i < k, let Hi = Gv \ Gvi = w LAw ∪ SAw , where w 6= v is a node lying on the path from v to vi . Let ρi := (ρ ∩ 2vi ) \ U(Hi ),
(3)
which is a rectangle obtained by clipping ρ from top and bottom (e.g. hatched rectangle in Figure 9 is ρ2 ). By construction, ρi is compatible with the right virtual wall of 2vi . Using (1) and (3), µ(Gv , ρ) =
k X i=0
=
k X i=0
µ(Hi , ρ ∩ 2vi ) + µ(Gvi , Hi , ρ ∩ 2vi ) µ(Hi , ρ ∩ 2vi ) + µ(Gvi , ρi ) .
⊖
The first term can be computed in O(log n) time. As for the second term, if vi is a leaf of T, then Gvi = ∅ and µ(Gvi , ρ) = 0, so assume that vi is an interior node of T. If vi ∈ D(v), then vi is ⊖type and ρi is compatible with L E+ vi (e.g. v1 , v3 in Figure 9). We can therefore compute µ(Gvi , ρi ) in O(log n) time, as described above. If vi 6∈ D(v), i.e., vi is -type and 1 ≤ i < k (e.g. v2 in Figure 9), then it can be proved that Bvi ∪ Fvi = ∅ and + µ(Gvi , ρi ) = µ(Cvi , ρi ) = µ(Cvi , P+ vi , ρi ) = C A (v, ρi ).
By Lemma 5.1, µ(Gvi , ρi ) can be computed in O(log n) time. Putting everything together, µ(Gv , ρ) can be computed in O(log2 n) time, and A DJUST WALL thus takes O(log2 n) time. L EMMA 4.1. A DJUST WALL and T IGHTEN W INDOW procedures take O(log2 n) time, and P USH D OWN takes O(1) time.
CORNER DATA STRUCTURE
Let 2 := [0, 0] × [a, b] be a rectangle with p0 , . . . , p3 as its vertices. Without loss of generality, assume that a ≥ b. Let C be a set of corner squares with respect to 2, let B be a set of rectangles, called bars, whose horizontal edges lie on the top and bottom edges of 2. Let I be the set of x-projections of bars in B. We describe a data structure that supports the following operations: insertion or deletion of a corner or a bar; for a query rectangle ρ ⊆ 2, compute C A(ρ) = µ(C, B, ρ). We assume that we have the wall data structure, described in Section 2, so that for any interval I, ε(I) = |I \ U(I)| can be computed in O(log n) time. For a corner square S, let q(S) be the vertex of S that lies inside 2. For 0 ≤ i ≤ 3, let Ci ⊆ C be the set of corners that contain the vertex pi of 2, and let Qi = {q(S) | S ∈ Ci }; see Figure 10 (ii). We call a (rectilinear) polygon staircase if it consists of a rectilinear chain that is both x- and y-monotone (in each coordinate it can be either increasing or decreasing), and the endpoints of the chain are connected together by a horizontal and a vertical edge; see Figure 10 (i). The common endpoint of the horizontal and the vertical edge is called the apex of the polygon. For each i, U(Ci , 2) is a staircase polygon whose apex is pi . We construct the following data structures: • We preprocess C0 and B into a data structure Ψ0 so that, for a query rectangle ρ ⊆ 2, µ(C0 , B, ρ) can be computed in O(log n) time. We construct a similar data structure Ψ1 on C1 and B for computing µ(C1 , B, ρ). • We preprocess C0 , C2 , B into a data structure Ψ2 so that, for a rectangle ρ ⊆ 2, µ(C2 , B ∪ C0 , ρ) can be computed in O(log n) time. We construct a similar data structure Ψ3 on C1 , C3 , B for computing µ(C3 , B ∪ C1 , ρ). Before describing Ψ0 , . . . , Ψ3 , we explain how they can be used to compute C A(ρ) in O(log n) time. For 0 ≤ i ≤ 3, let Γi be the rectilinear chain of U(Ci ); Figure 10 (i). Γi and Γi+1 intersect in at most one point, but Γi and Γi+2 may intersect several times. For each i, we define a point ri , as follows. If Γi and Γi+1 (mod 4) do not intersect, then ri is the endpoint of Γi that lies on the edge pi pi+1 , otherwise ri is the (orthogonal) projection of the intersection point of Γi and Γi+1 on the edge pi pi+1 . Finally, let 2i be the rectangle formed by ri−1 and ri as two of its diagonally opposite vertices; see Figure 10 (iv) and (v). The interiors of 2i and 2i+1 are disjoint but those of 2i and 2i+2 may intersect. However, if 20 and 22 intersect, then 21 and 23 are disjoint, and vice-versa; Figure 10 (iv). Note that 20 , . S . . , 23 may not cover 2; Figure 10 (v). Nevertheless, U(C, 2) ⊆ 3i=0 2i and Γi appears on ∂U(C, 2) only inside 2i ; see Figure S 10. The following equality is straightforward: U(C, 2) = 3i=0 U(Ci , 2i ). For a rectangle ρ, we set ρi = ρ ∩ 2i . Then U(C, B, ρ) =
3 [
U(Ci , B, ρi ).
i=0
Since the interiors of 2i and 2i+1 are disjoint, µ(C, B, ρ)
=
µ(C0 , B, ρ0 ) + µ(C1 , B, ρ1 ) + µ(C2 , C0 ∪ B, ρ2 ) + µ(C3 , C1 ∪ B, ρ3 ) .
Let Vi ⊆ Qi be the subset of points that appear as (convex) vertices of U(Ci ). We maintain Vi sorted from left to right. By performing a binary search on Vi and Vi+1 , ri can be computed in O(log n) time. We can then compute ρ0 , . . . , ρ3 in additional O(1) time. Each of the terms on the right hand side of the above equation can be computed in O(log n) time using Ψ0 , . . . , Ψ3 , therefore
p3
p2
Γ0 p1
(i)
(ii)
r2
23
Γ2
Γ3
p0
r2
p2
p3
22 r1
r1
Γ1
p0
23
22
r3 20
21
p1
r3
20
r0
(iii)
21 r0
(iv)
(v)
Figure 10. (i) Staircase polygons. (ii) Four different types of corner squares. (iii) rectilinear chains Γ0 , . . . , Γ3 . (iv) 20 , . . . , 23 cover 2; 20 and 22 intersect. (v) 20 , . . . , 23 that do not cover 2.
µ(C, B, ρ) can also be computed in O(log n) time. We now describe Ψ0 and Ψ2 ; Ψ1 and Ψ3 are symmetric to Ψ0 and Ψ2 , respectively. Although Ψ0 is a special case of Ψ2 , we first describe Ψ0 and then describe how we extend it to construct Ψ2 . Preprocessing C0 and B. For simplicity, we replace each square C in C0 into a quadrant with apex q(C) and containing p0 ; Let Q be the resulting set of quadrants. We describe the data structure to preprocess Q and B to compute µ(Q, B, ρ) for a query rectangle ρ ⊆ 2. U(Q) can be maintained using the data structure by Overmars and van Leeuwen [14] but the interaction of Q and B makes the maintenance of U(Q, B) difficult: a bar in B may intersect many edges of ∂U(Q), and many bars in B may intersect the same edge of U(Q); see Figure 11. The insertion or deletion of a corner square or a bar may thus cause many changes in the boundary structure of U(Q, B), making it expensive to update U(Q, B). We circumvent this problem by maintaining the edges of U(Q, B) implicitly, as described below. B2
B1 a b
c
B4
d e f g
B3
Figure 11. Interaction of Q and B. Γ = ha, . . . , gi; b, c, d do not appear ˆ = hˆ on Σ(Γ, B); Σ a, eˆ, fˆ.ˆ g i.
Operations on a monotone chain. For a rectilinear xy-monotone chain Γ ⊂ 2, let P(Γ) be the staircase polygon formed by the set of points in 2 that lie below Γ. We represent Γ by its sequence of horizontal edges, and we will use Γ to denote both the chain as well as the sequence of its horizontal edges. For a subset B ⊆ B, let Σ(Γ, B) be the sequence of horizontal edges of Γ \ U(B). An edge of Γ may get split into many edges in Σ(Γ, B). We represent Σ(Γ, B) implicitly: For an edge e ∈ Γ, define a 4tuple eˆ = hx− (e), x+ (e), y(e), ε(e)i, where x− (e), x+ (e) are the x-coordinates of the left and right endpoints of e, y(e) is the ycoordinate of e, and ε(e) = |e \ U(B)| is the exposed portion of e on Σ; eˆ is the representation of e. We say that e appears in Σ ˆ = hˆ if ε(e) > 0. We represent Σ as the sequence Σ e1 , eˆ2 , . . . , i, where e1 , e2 , . . . is the subsequence of P the edges of Γ that appear in Σ. Note that Area(P(Γ) \ U(B)) = eˆ∈Σ ε(e) · y(e). Abusing the notation slightly, we do not distinguish between Σ and its representation. For A ⊂ R2 , let X(A) denote the x-projection of A. We perform the following operations on Σ. (i) S PLIT(Σ, a): Given a rectilinear chain Σ and a value a ≥ 0, let Σ+ = he ∈ Σ | y(e) > ai and Σ− = he ∈ Σ |
y(e) < ai. Split Σ into Σ− , Σ+ . If Σ contains an edge e with y(e) = a, then the procedure returns eˆ separately. See Figure 12. (ii) M ERGE(Σ− , Σ+ , e): Let Σ− , Σ+ be two xy-monotone rectilinear and e a horizontal edge so that Σ = Σ− ◦ hei ◦ Σ+ is also a xy-monotone rectilinear chain. Given Σ− , Σ+ , and eˆ, compute Σ. If the edge e is not specified or ε(e) = 0, the procedure computes Σ1 ◦ Σ2 . (iii) Q UERY(Σ, I): Given a chain Σ and an interval I, compute X α(Σ, I) = ε(e) · y(e). e∈Σ,X(e)⊆I
Suppose ei , . . . , ej are the edges of Σ whose x-projections are contained in I. Then α(Σ, I) is the area of (P(Γ) \ U(B)) ∩ RI , where RI = [x− (ei ), x+ (ej )] × [0, b].
11111 00000 00000 11111 00000 11111 00000 11111 00000 11111 00000 11111 00000 11111 00000 11111 00000 11111 00000 11111 00000 11111 00000 11111 00000 11111 00000 11111
RI
I
11 00 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 00 11 11 00 11
S PLIT a
M ERGE
1111 0000 0000 1111 0000 1111 0000 1111 0000 1111 0000 1111 0000 1111 0000 1111 0000 1111 0000 1111 0000 1111 0000 1111 0000 1111
000 111111 000 000 111 000 111 000 111 000 111 000 111 000 111 000 111 000 111 000 111 000 111 000 111 000 e111 000 111 000 111 000 111 000 111 000 111 Σ− 000 111 000 111 000 111 000 000111 111
Σ+111 000 000111
a
Figure 12. Operations on chains; I is a query interval, darker region is (P(Γ) \ U(B)) ∩ RI (left).
By maintaining Σ in a height-balanced tree, which we denote by H(Σ), each of these operations can be performed in O(log n) time. Using the wall data structure, for an edge e of Γ, ε(e) can be computed in O(log n) time. δu
ψ(u, p)
111 1 000 0 000 111 0 1 000 111 0 ℓu 1 ⊓ 000 111 0 1 pu 000 111 0 1 000 111 0 000 1 111 0 1 000 1 111 0 br(u) 000 111 0 1 000 111 0 1 ˜ 1 000 111 Σ(w) 000 0 111 0 1 000 1 111 0 δw
11 11 00 00 00 11 00 00 11 11 00 11 00 11 00 11 00 11 00 11 p 00 11 00 00 11 11 00 11 00 11 11 00 00 11 00 11 00 11 00 11 00 11 00 00 11 11 00 11 00 11 11 00
δz
Figure 13. A node u and its children w and z. B∗ (u) is the hatched area and U(Qu , B∗u , 2u ) is the darker region.
Structure of Ψ0 . Let Q = {q(S) | S ∈ Q}, and let X ⊆ [0, a] be the set of points in R that are x-coordinates of points in Q and vertices in B. X partitions the interval (0, a] into a set of atomic intervals; we regard each of these intervals semiclosed — the right
endpoint lies in the interval but the left endpoint does not lie in the interval. We build a height-balanced binary tree T on these atomic intervals. Each node v ∈ T is associated with an interval + δv = (x− v , xv ] and a rectangle 2v = δv × [0, b]. If v is the i-th leaf of T, then δv is the i-th (semiclosed) atomic interval. For an interior node v, let δv = δL(v) ∪ δR(v) , 2v = δv × [0, b], and let ξv be the common endpoint of δL(v) and δR(v) . Let ℓv : x = ξv be the separator line at v. For an interval I ⊆ [0, a], let Ξ(I) be the set of nodes u ∈ T for which clδu ⊆ I and clδp(u) 6⊆ I, and let N(I) be the set of ancestors of nodes in Ξ(I). Similarly, we define Ξ(B) and N(B) for a bar in B. For a node u, let Qu = 2u ∩ Q, Qu = {S ∈ Q | q(S) ∈ Qu }, Bu ⊆ B (resp. B∗u ⊆ B) be the set of bars B such that u ∈ Ξ(B) (resp. u ∈ Ξ(B) ∪ N(B)), and p⊓ u be the point of Qu with the maximum y-coordinate. Let β(u) = |X(B∗ (u))∩δv |, which is |δu | if |Bu | > 0, and β(L(u))+β(R(u)) otherwise. U(Qu , 2u ) is a staircase polygon. Let Γ(u) be its rectilinear chain. Set Σ(u) = Σ(Γ(u), B∗u ), the portion of ∂U(Qu , B∗u , 2u ) that lies on Γ(u). For an interior node u, let br(u) be the edge ˜ of Γ(u) that intersects ℓu , and let Σ(u) = Σ(u) \ Σ(p(u)) (Figure 15 (i)). If ε(br(u)) > 0, then br(u) ∈ Σ(u). For the root rt of ˜ ˜ T, Σ(rt) = Σ(rt). At each node u we store: Σ(u), br(u), β(u), |Bu |, and p⊓ (u). For a node u and for a point p, let ψ(u, p) denote the rightmost point in Qu whose y-coordinate is greater than that of p; ψ(u, p) can be computed in O(log n) time. Updating Ψ0 . Using the following two procedures, each of which takes O(log n) time, a corner or a bar can be inserted/deleted in O(log2 n) time. Let u be a node of T, and let w and z be the left and right child of u, respectively. ˜ ˜ A SCEND(u): It computes Σ(u), Σ(w), Σ(z) from Σ(w) and Σ(z). ˜ If β(u) = δu , then we set Σ(u) = ∅, Σ(w) = Σ(w), and ˜ Σ(z) = Σ(z), and we stop. So assume that β(u) < δu . Let pR = p⊓ (z). We compute pL = ψ(w, pR ). If pL does not exist, then the left endpoint of br(u) is (x− u , y(pR )), otherwise it is (x(pL ), y(pR )). The right endpoint of br(u) is pR . Let I = [x(pL ), x(pR )] be the projection of br(u). We compute ε(br(u)), in O(log n) time. If ε(br(u)) = 0, we set br(u) = ∅. Next, we call S PLIT(Σ(w), y(pR )). Let Σ− , Σ+ be the subchains returned by the procedure. We set ˜ Σ(w) = Σ+ . If Σ(z) contains the horizontal edge at height y(pR ), we delete it, and we set Σ(u) = M ERGE(Σ− , Σ(z), br(u)).
111 000 100 2w1 000 111 000 111 1010 000 111 p⊓ w 000 111 100 000 111 p1 000 111 1010L 000 111 000 111 1010 000 111 000 111 10 000 111 000 1010 111
11 00 00 2z 11 00 11 00 ℓu 11 00 11 00 11 ⊓ p11 00 z 00 11 00 11 00 11 00 11 00 11 00 11
11 00 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11
2u
111 000 1010 000 111 000 111 1010 000 111 p⊓ u 000 111 1010 000 111 000 111 1010 br(u) 000 111 000 111 100 ˜ 1 000 111 Σ(w) 000 111 10 000 111 000 1010 111
11 00 00 11 00 11 00 11 00 11 00 11 00 p11 R 00 11 00 11 00 11 00 11 00 11 00 11
11 00 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11
e−
(x+ , y + ) eH
Iy−
00 11 ρ 00 11 00 11 000000 00111111 11 00111111 11 R− 11 000000 00 eL 00 11 000000 111111 00 11 P 00000000 11 00111111 11
00 11111111 000000 00 000000 111111 R 11 000000 111111
11 00 00 11
(x− , y − )
Ix−
I
Ix+
R+
Iy+ e+
Figure 15. Answering a query on Ψ0 .
Answering a query. Let ρ = [x− , x+ ] × [y − , y + ] ⊆ 2 be a query rectangle. We wish to compute µ(Q, B, ρ). Let rt be the root of T, and let Γ = Γ(rt) and Σ = Σ(rt). By searching with ρ in H(Σ), stored at rt, we find the highest (resp. lowest) edge eH (resp. eL ) of Γ that appears in Σ and that lies completely inside ρ. Let P be the staircase polygon formed by the edges eH , . . . , eL of Σ. Let e− be the predecessor of eH and e+ the successor of eL in Σ. Set I = [x− (eH ), x+ (eL )], R = I ×[0, y − ], Ix− = X(e− )∩[x− , x+ ], Ix+ = X(e+ ) ∩ [x− , x+ ], Iy− = [0, y(e− )] ∩ [y − , y + ], and Iy+ = [0, y(e+ )] ∩ [y − , y + ]. Set R− = Ix− × Iy− and R+ = I + × Iy+ . See Figure 15 (ii). Let A (resp. A− , A+ ) denote the area of R (resp. R− , R+ ) not lyig in U(B). Then µ(Q, B, ρ)
= α(Σ, I) − A + A− + A+ = α(Σ, I) − ε(I)yL + ε(Ix− )|Iy− | + ε(I + )|Iy+ |.
Each of these quantities can be computed in O(log n) time either querying H(Σ) or the wall data structure. Hence, µ(Q, B, ρ) can be computed in O(log n) time. Preprocessing C0 , C2 , and B. The overall structure of Ψ2 is similar to that of Ψ0 , but it is more complex because of the interaction between C0 and C2 . To keep the presentation analogous to Ψ0 , we describe the data structure for computing µ(C0 , C2 ∪ B, ρ); µ(C2 , C0 ∪ B, ρ) can be computed by reversing the directions of (+x)- and (+y)-axes and constructing the data structure in this new coordinate frame. Let Q+ (resp. Q− ) be the set of quadrants corresponding to the squares in C0 and C2 . For a query rectangle ρ ⊆ 2, we wish to compute µ(Q+ , Q− ∪ B, ρ). We begin by describing the representation of U(Q+ , Q− ∪ B, 2) and the operations performed on it.
(i)
(ii)
000 111 111 000 000 111 000 111 000 111 000 111 000 111 000 111 000 111 000 111 000 111 000 111 000 111 000 111 000 111 000 111 000 111 000 111
(iii)
00 11 11 00 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11 00 11
Figure 16. Representation of Σ(Q+ , Q− , B). (i) U(Q+ ). (ii) U(Q+ , Q− ), black (resp. white) points are vertices of squares in Q+ (resp. Q− ), and gray vertices are mixed vertices. (iii) U(Q+ , Q− ∪ B); edges inside the ellipse form a superedge.
Figure 14. A SCEND(u) and D ESCEND(u) procedures.
D ESCEND(u): Assuming Σ(u) is at our disposal, this procedure computes Σ(w) and Σ(z). If β(u) = δu , then we set Σw = ˜ w and Σ(z) = Σ(z) ˜ Σ and stop. Otherwise, let pR = p⊓ (z). We call S PLIT(Σ(u), y(pR )). Let Σ− , Σ+ be the two chains returned by this procedure; if ε(br(u)) > 0, then the proce˜ dure also returns br(u). We set Σ(w) = M ERGE(Σ− , Σ(w)). Let e = br(u) ∩ 2z . We compute ε(e) and set Σ(z) = M ERGE(∅, Σ+ , e).
Representation of U(Q+ , Q− ∪ B, 2). U(Q+ , Q− , 2) is a staircase polygon, with four types of vertices: vertices of Q+ , vertices of Q− , intersection points of two edges of Q+ (or Q− ), and intersection points of an edge of ∂U(Q+ ) with an edge of ∂U(Q− ). The first and the last type are convex vertices, the second type are reflex vertices, and the third type can be both convex and reflex. We refer to the last ones as mixed vertices. An edge of U(Q+ , Q− ) is called purple if one of its endpoints is a mixed vertex, blue if it is an edge of U(Q+ ), and red if it is an edge of U(Q− ). Let
Γ(Q+ , Q− ) be the rectilinear chain of U(Q+ , Q− , 2). As in the previous case, Γ is represented as the sequence of its horizontal edges. Let Σ := Σ(Q+ , Q− , B) be the sequence of edges of Γ(Q+ , Q− ) \ U(B), and let Σ(Q− , B) be the sequence of edges in Γ(Q− ) \ U(B). Σ is represented implicitly using Σ(Q− , B): Let E = he1 , e2 , . . . , ek i be the sequence of edges of Γ(Q+ , Q− ) that appear in Σ. Every red edge in E also appears in Σ(Q− , B). We call a contiguous subsequence ei , ei+1 , . . . , ej of red edges in E as a super-edge σij (see Figure 16 (iii)), and represent it as a five-tuple: σ ˆij = hx− (σij ), x+ (σij ), y − (σij ), y + (σij ), ε(σij )i,
where x− (σij ) = x− (ei ), x+ (σij ) = x+ (ej ), y − (σij ) = y(ej ), Pj y + (σij ) = y(ei ), and ε(σij ) = l=i ε(el )y(el ). We represent ˆ = hγ1 , γ2 , . . . , γs i as a sequence of blue edges, purple Σ = Σ edges, and super-edges; either γi = eˆj for a blue or purple edge in E or γi = σ ˆjl for a super-edge σjl . By definition, a super-edge may consist of a single red edge. We assume that each super-edge ˆ is maximal in the sense that if γi is a super-edge then γi+1 is in Σ not a super-edge. By storing both Σ(Q+ , Q− , B) and Σ(Q− , B) in height balanced binary trees, each of the following operations can be performed in O(log n) time. (i) S PLIT(Σ, a): Given a rectilinear chain Σ and a value a ≥ 0, let Σ+ = he ∈ Σ | y(e) > ai and Σ− = he ∈ Σ | y(e) < ai. Split Σ into Σ− , Σ+ . If Σ contains an edge e with y(e) = a, then the procedure returns eˆ separately. If ˆ with a ∈ [y − (σjl ), y + (σjl )], there is a super-edge σjl in Σ then the procedure splits the super-edge into two: one of which lies above a and appears at the end of Σ− , and the other lies below a and appears at the beginning of Σ+ . (ii) M ERGE(Σ− , Σ+ , e): Let Σ− , Σ+ be two xy-monotone rectilinear chains and e a horizontal edge so that Σ = Σ− ◦hei◦ Σ+ is also a xy-monotone rectilinear chain. Given Σ− , Σ+ , and eˆ, compute Σ. If the edge e is not specified or ε(e) = 0, the procedure computes Σ− ◦ Σ+ . Finally, if there are two ˆ then we merge them. consecutive super-edges in Σ, (iii) Q UERY(Σ, I): Given a chain Σ and an interval I, compute X α(Σ, I) = ε(e) · y(e). e∈Σ,X(e)⊆I
Since some of the edges of Σ are implicitly represented as super-edges, computing α(Σ, I) now is more involved than for Ψ0 but nevertheless can be computed in O(log n) time. Structure of Ψ2 . Let Q− = {q(S) | S ∈ Q− } and Q+ = {q(S) | S ∈ Q+ }. Let X ⊆ [0, a] be the set of x-coordinates of points in Q+ ∪ Q− and edges in B. We construct a binary tree T as − ∗ for Ψ0 . For each node u, we define Q+ u , Qu , Bu , Bu in a man− − + + ner similar to Ψ0 . Let Γ (u) = Γ(Qu , Qu ), Γ (u) = Γ(Q− u ), − ∗ − − ∗ Σ+ (u) = Σ(Q+ u , Qu , Bu ), and Σ (u) = Σ(Qu , Bu ). In ad˜ + (u) = Σ+ (u) \ Σ+ (p(u)) and Σ ˜ − (u) = Σ− (u) \ dition, let Σ − + − Σ (p(u)). We define br (u) (resp. br (u)) be the edge of Γ+ (u) (resp. Γ− (u)) that intersects ℓu . At each node u, we store compact ˜ + (u) and Σ ˜ − (u), as described above, and we representations of Σ + − also store br (u) and br (u). As for Ψ0 , we again describe A SCEND, D ESCEND procedures to compute Σ+ (u) from Σ+ (w) and Σ+ (z), and vice-versa, where w and z are the children of v. Each of them takes O(log n) time, and they are used to insert and delete corners and bars. Putting everything together, we conclude the following:
L EMMA 5.1. A set C of corner squares and a set B of bars can be maintained in a dynamic data structure under insertion/deletion so that for a query rectangle ρ, µ(C, B, ρ) can be computed in O(log n) time. The insertion and deletion operations take O(log2 n) time.
References [A1] http://en.wikipedia.org/wiki/Klee’s measure problem [A2] P. K. Agarwal, H. Kaplan, and M. Sharir, Computing the volume of the union of cubes, Proc. 23rd Annu. Sympos. Comput. Geom., 2007, 294–301. [A3] M. de Berg, M. van Kreveld, M. Overmars and O. Schwarzkopf, Computational Geometry: Algorithms and Applications, 2nd edition, Springer Verlag, Heidelberg, 2000. [A4] J. L. Bentley, Algorithms for Klee’s rectangle problems. Unpublished notes, Computer Science Department, Carnegie Mellon University, 1977. [A5] J.D. Boissonnat, M. Sharir, B. Tagansky and M. Yvinec, Voronoi diagrams in higher dimensions under certain polyhedral distance functions, Discrete Comput. Geom. 19 (1998), 485– 519. [A6] K. Bringmann, Klee’s measure problem on fat boxes, Proc. 26th Annu. Sympos. Comput. Geom., 2010, to appear. [A7] T. M. Chan, A (slightly) faster algorithm for Klee’s measure problem, Comput. Geom.: Theory Appls. 43 (2010), 243-250. [A8] E. Chen and T. M. Chan, Space-efficient algorithms for Klee’s measure problem, Proc. 17th Canadian Conf. Comput. Geom., 2005. [A9] B. S. Chlebus, On the Klee’s measure problem in small dimensions, Proc. 25th Conf. Current Trends in Theory and Practice of Informatics, 1998, 304-311. [A10] M. Dickerson, C. A. Duncan, and M. T. Goodrich, K-d trees are better when cut on the longest side, Proc. 8th European Sympos. Algorithms, 2000, 179190. [A11] M. L. Fredman S and B. Weide, The complexity of computing the measure of [ai , bi ], Commun. ACM 21 (1978), 540–544. [A12] J. van Leeuwen and D. Wood, The measure problem for rectangular ranges in d-space, J. Algorithms 2 (1981), 282–300. S [A13] V. Klee, Can the measure of [ai , bi ] be computed in less than O(n log n) steps? Amer. Math. Monthly 84 (1977), 284– 285. [A14] M. Overmars and J. Leeuwen, Maintenance of configurations in the plane, J. Comput. Syst. Sci. 23 (1981), 166–204. [A15] M. Overmars and C.K. Yap, New upper bounds in Klee’s measure problem, SIAM J. Comput. 20 (1991), 1034–1045. [A16] S. Suzuki and T. Ibaki, An average running time analysis of a backtracking algorithm to calculate the measure of the union of hyperrectangles in d dimensions, Proc. 16th Canadian Conf. Comput. Geom., 2004, 196–199.