1// Copyright 2005 Google Inc. All Rights Reserved.
2
3#ifndef UTIL_GEOMETRY_S2CELLID_H_
4#define UTIL_GEOMETRY_S2CELLID_H_
5
6#include <iostream>
7using std::ostream;
8using std::cout;
9using std::endl;
10
11#include <string>
12using std::string;
13
14#include <vector>
15using std::vector;
16
17#include "base/basictypes.h"
18#include "base/logging.h"
19#include "base/port.h" // for HASH_NAMESPACE_DECLARATION_START
20#include "s2.h"
21#include "util/math/vector2.h"
22
23class S2LatLng;
24
25// An S2CellId is a 64-bit unsigned integer that uniquely identifies a
26// cell in the S2 cell decomposition. It has the following format:
27//
28// id = [face][face_pos]
29//
30// face: a 3-bit number (range 0..5) encoding the cube face.
31//
32// face_pos: a 61-bit number encoding the position of the center of this
33// cell along the Hilbert curve over this face (see the Wiki
34// pages for details).
35//
36// Sequentially increasing cell ids follow a continuous space-filling curve
37// over the entire sphere. They have the following properties:
38//
39// - The id of a cell at level k consists of a 3-bit face number followed
40// by k bit pairs that recursively select one of the four children of
41// each cell. The next bit is always 1, and all other bits are 0.
42// Therefore, the level of a cell is determined by the position of its
43// lowest-numbered bit that is turned on (for a cell at level k, this
44// position is 2 * (kMaxLevel - k).)
45//
46// - The id of a parent cell is at the midpoint of the range of ids spanned
47// by its children (or by its descendants at any level).
48//
49// Leaf cells are often used to represent points on the unit sphere, and
50// this class provides methods for converting directly between these two
51// representations. For cells that represent 2D regions rather than
52// discrete point, it is better to use the S2Cell class.
53//
54// This class is intended to be copied by value as desired. It uses
55// the default copy constructor and assignment operator.
56class S2CellId {
57 public:
58 // Although only 60 bits are needed to represent the index of a leaf
59 // cell, we need an extra bit in order to represent the position of
60 // the center of the leaf cell along the Hilbert curve.
61 static int const kFaceBits = 3;
62 static int const kNumFaces = 6;
63 static int const kMaxLevel = S2::kMaxCellLevel; // Valid levels: 0..kMaxLevel
64 static int const kPosBits = 2 * kMaxLevel + 1;
65 static int const kMaxSize = 1 << kMaxLevel;
66
67 inline explicit S2CellId(uint64 id) : id_(id) {}
68
69 // The default constructor returns an invalid cell id.
70 inline S2CellId() : id_(0) {}
71 inline static S2CellId None() { return S2CellId(); }
72
73 // Returns an invalid cell id guaranteed to be larger than any
74 // valid cell id. Useful for creating indexes.
75 inline static S2CellId Sentinel() { return S2CellId(~uint64(0)); }
76
77 // Return a cell given its face (range 0..5), 61-bit Hilbert curve position
78 // within that face, and level (range 0..kMaxLevel). The given position
79 // will be modified to correspond to the Hilbert curve position at the
80 // center of the returned cell. This is a static function rather than a
81 // constructor in order to give names to the arguments.
82 static S2CellId FromFacePosLevel(int face, uint64 pos, int level);
83
84 // Return the leaf cell containing the given point (a direction
85 // vector, not necessarily unit length).
86 static S2CellId FromPoint(S2Point const& p);
87
88 // Return the leaf cell containing the given normalized S2LatLng.
89 static S2CellId FromLatLng(S2LatLng const& ll);
90
91 // Return the direction vector corresponding to the center of the given
92 // cell. The vector returned by ToPointRaw is not necessarily unit length.
93 S2Point ToPoint() const { return ToPointRaw().Normalize(); }
94 S2Point ToPointRaw() const;
95
96 // Return the center of the cell in (s,t) coordinates (see s2.h).
97 Vector2_d GetCenterST() const;
98
99 // Return the center of the cell in (u,v) coordinates (see s2.h). Note that
100 // the center of the cell is defined as the point at which it is recursively
101 // subdivided into four children; in general, it is not at the midpoint of
102 // the (u,v) rectangle covered by the cell.
103 Vector2_d GetCenterUV() const;
104
105 // Return the S2LatLng corresponding to the center of the given cell.
106 S2LatLng ToLatLng() const;
107
108 // The 64-bit unique identifier for this cell.
109 inline uint64 id() const { return id_; }
110
111 // Return true if id() represents a valid cell.
112 inline bool is_valid() const;
113
114 // Which cube face this cell belongs to, in the range 0..5.
115 inline int face() const;
116
117 // The position of the cell center along the Hilbert curve over this face,
118 // in the range 0..(2**kPosBits-1).
119 inline uint64 pos() const;
120
121 // Return the subdivision level of the cell (range 0..kMaxLevel).
122 int level() const;
123
124 // Return the edge length of this cell in (i,j)-space.
125 inline int GetSizeIJ() const;
126
127 // Return the edge length of this cell in (s,t)-space.
128 inline double GetSizeST() const;
129
130 // Like the above, but return the size of cells at the given level.
131 inline static int GetSizeIJ(int level);
132 inline static double GetSizeST(int level);
133
134 // Return true if this is a leaf cell (more efficient than checking
135 // whether level() == kMaxLevel).
136 inline bool is_leaf() const;
137
138 // Return true if this is a top-level face cell (more efficient than
139 // checking whether level() == 0).
140 inline bool is_face() const;
141
142 // Return the child position (0..3) of this cell's ancestor at the given
143 // level, relative to its parent. The argument should be in the range
144 // 1..kMaxLevel. For example, child_position(1) returns the position of
145 // this cell's level-1 ancestor within its top-level face cell.
146 inline int child_position(int level) const;
147
148 // Methods that return the range of cell ids that are contained
149 // within this cell (including itself). The range is *inclusive*
150 // (i.e. test using >= and <=) and the return values of both
151 // methods are valid leaf cell ids.
152 //
153 // These methods should not be used for iteration. If you want to
154 // iterate through all the leaf cells, call child_begin(kMaxLevel) and
155 // child_end(kMaxLevel) instead.
156 //
157 // It would in fact be error-prone to define a range_end() method,
158 // because (range_max().id() + 1) is not always a valid cell id, and the
159 // iterator would need to be tested using "<" rather that the usual "!=".
160 inline S2CellId range_min() const;
161 inline S2CellId range_max() const;
162
163 // Return true if the given cell is contained within this one.
164 inline bool contains(S2CellId const& other) const;
165
166 // Return true if the given cell intersects this one.
167 inline bool intersects(S2CellId const& other) const;
168
169 // Return the cell at the previous level or at the given level (which must
170 // be less than or equal to the current level).
171 inline S2CellId parent() const;
172 inline S2CellId parent(int level) const;
173
174 // Return the immediate child of this cell at the given traversal order
175 // position (in the range 0 to 3). This cell must not be a leaf cell.
176 inline S2CellId child(int position) const;
177
178 // Iterator-style methods for traversing the immediate children of a cell or
179 // all of the children at a given level (greater than or equal to the current
180 // level). Note that the end value is exclusive, just like standard STL
181 // iterators, and may not even be a valid cell id. You should iterate using
182 // code like this:
183 //
184 // for(S2CellId c = id.child_begin(); c != id.child_end(); c = c.next())
185 // ...
186 //
187 // The convention for advancing the iterator is "c = c.next()" rather
188 // than "++c" to avoid possible confusion with incrementing the
189 // underlying 64-bit cell id.
190 inline S2CellId child_begin() const;
191 inline S2CellId child_begin(int level) const;
192 inline S2CellId child_end() const;
193 inline S2CellId child_end(int level) const;
194
195 // Return the next/previous cell at the same level along the Hilbert curve.
196 // Works correctly when advancing from one face to the next, but
197 // does *not* wrap around from the last face to the first or vice versa.
198 inline S2CellId next() const;
199 inline S2CellId prev() const;
200
201 // This method advances or retreats the indicated number of steps along the
202 // Hilbert curve at the current level, and returns the new position. The
203 // position is never advanced past End() or before Begin().
204 S2CellId advance(int64 steps) const;
205
206 // Like next() and prev(), but these methods wrap around from the last face
207 // to the first and vice versa. They should *not* be used for iteration in
208 // conjunction with child_begin(), child_end(), Begin(), or End(). The
209 // input must be a valid cell id.
210 inline S2CellId next_wrap() const;
211 inline S2CellId prev_wrap() const;
212
213 // This method advances or retreats the indicated number of steps along the
214 // Hilbert curve at the current level, and returns the new position. The
215 // position wraps between the first and last faces as necessary. The input
216 // must be a valid cell id.
217 S2CellId advance_wrap(int64 steps) const;
218
219 // Iterator-style methods for traversing all the cells along the Hilbert
220 // curve at a given level (across all 6 faces of the cube). Note that the
221 // end value is exclusive (just like standard STL iterators), and is not a
222 // valid cell id.
223 inline static S2CellId Begin(int level);
224 inline static S2CellId End(int level);
225
226 // Methods to encode and decode cell ids to compact text strings suitable
227 // for display or indexing. Cells at lower levels (i.e. larger cells) are
228 // encoded into fewer characters. The maximum token length is 16.
229 //
230 // ToToken() returns a string by value for convenience; the compiler
231 // does this without intermediate copying in most cases.
232 //
233 // These methods guarantee that FromToken(ToToken(x)) == x even when
234 // "x" is an invalid cell id. All tokens are alphanumeric strings.
235 // FromToken() returns S2CellId::None() for malformed inputs.
236 string ToToken() const;
237 static S2CellId FromToken(string const& token);
238
239 // Creates a debug human readable string. Used for << and available for direct
240 // usage as well.
241 string ToString() const;
242
243 // Return the four cells that are adjacent across the cell's four edges.
244 // Neighbors are returned in the order defined by S2Cell::GetEdge. All
245 // neighbors are guaranteed to be distinct.
246 void GetEdgeNeighbors(S2CellId neighbors[4]) const;
247
248 // Return the neighbors of closest vertex to this cell at the given level,
249 // by appending them to "output". Normally there are four neighbors, but
250 // the closest vertex may only have three neighbors if it is one of the 8
251 // cube vertices.
252 //
253 // Requires: level < this->level(), so that we can determine which vertex is
254 // closest (in particular, level == kMaxLevel is not allowed).
255 void AppendVertexNeighbors(int level, vector<S2CellId>* output) const;
256
257 // Append all neighbors of this cell at the given level to "output". Two
258 // cells X and Y are neighbors if their boundaries intersect but their
259 // interiors do not. In particular, two cells that intersect at a single
260 // point are neighbors.
261 //
262 // Requires: nbr_level >= this->level(). Note that for cells adjacent to a
263 // face vertex, the same neighbor may be appended more than once.
264 void AppendAllNeighbors(int nbr_level, vector<S2CellId>* output) const;
265
266 /////////////////////////////////////////////////////////////////////
267 // Low-level methods.
268
269 // Return a leaf cell given its cube face (range 0..5) and
270 // i- and j-coordinates (see s2.h).
271 static S2CellId FromFaceIJ(int face, int i, int j);
272
273 // Return the (face, i, j) coordinates for the leaf cell corresponding to
274 // this cell id. Since cells are represented by the Hilbert curve position
275 // at the center of the cell, the returned (i,j) for non-leaf cells will be
276 // a leaf cell adjacent to the cell center. If "orientation" is non-NULL,
277 // also return the Hilbert curve orientation for the current cell.
278 int ToFaceIJOrientation(int* pi, int* pj, int* orientation) const;
279
280 // Return the lowest-numbered bit that is on for this cell id, which is
281 // equal to (uint64(1) << (2 * (kMaxLevel - level))). So for example,
282 // a.lsb() <= b.lsb() if and only if a.level() >= b.level(), but the
283 // first test is more efficient.
284 uint64 lsb() const { return id_ & -id_; }
285
286 // Return the lowest-numbered bit that is on for cells at the given level.
287 inline static uint64 lsb_for_level(int level) {
288 return uint64(1) << (2 * (kMaxLevel - level));
289 }
290
291 private:
292 // This is the offset required to wrap around from the beginning of the
293 // Hilbert curve to the end or vice versa; see next_wrap() and prev_wrap().
294 static uint64 const kWrapOffset = uint64(kNumFaces) << kPosBits;
295
296 // Return the i- or j-index of the leaf cell containing the given
297 // s- or t-value. Values are clamped appropriately.
298 inline static int STtoIJ(double s);
299
300 // Return the (face, si, ti) coordinates of the center of the cell. Note
301 // that although (si,ti) coordinates span the range [0,2**31] in general,
302 // the cell center coordinates are always in the range [1,2**31-1] and
303 // therefore can be represented using a signed 32-bit integer.
304 inline int GetCenterSiTi(int* psi, int* pti) const;
305
306 // Given (i, j) coordinates that may be out of bounds, normalize them by
307 // returning the corresponding neighbor cell on an adjacent face.
308 static S2CellId FromFaceIJWrap(int face, int i, int j);
309
310 // Inline helper function that calls FromFaceIJ if "same_face" is true,
311 // or FromFaceIJWrap if "same_face" is false.
312 inline static S2CellId FromFaceIJSame(int face, int i, int j, bool same_face);
313
314 uint64 id_;
315} PACKED; // Necessary so that structures containing S2CellId's can be PACKED.
316DECLARE_POD(S2CellId);
317
318inline bool operator==(S2CellId const& x, S2CellId const& y) {
319 return x.id() == y.id();
320}
321
322inline bool operator!=(S2CellId const& x, S2CellId const& y) {
323 return x.id() != y.id();
324}
325
326inline bool operator<(S2CellId const& x, S2CellId const& y) {
327 return x.id() < y.id();
328}
329
330inline bool operator>(S2CellId const& x, S2CellId const& y) {
331 return x.id() > y.id();
332}
333
334inline bool operator<=(S2CellId const& x, S2CellId const& y) {
335 return x.id() <= y.id();
336}
337
338inline bool operator>=(S2CellId const& x, S2CellId const& y) {
339 return x.id() >= y.id();
340}
341
342inline bool S2CellId::is_valid() const {
343 return (face() < kNumFaces && (lsb() & GG_ULONGLONG(0x1555555555555555)));
344}
345
346inline int S2CellId::face() const {
347 return id_ >> kPosBits;
348}
349
350inline uint64 S2CellId::pos() const {
351 return id_ & (~uint64(0) >> kFaceBits);
352}
353
354inline int S2CellId::GetSizeIJ() const {
355 return GetSizeIJ(level());
356}
357
358inline double S2CellId::GetSizeST() const {
359 return GetSizeST(level());
360}
361
362inline int S2CellId::GetSizeIJ(int level) {
363 return 1 << (kMaxLevel - level);
364}
365
366inline double S2CellId::GetSizeST(int level) {
367 // Floating-point multiplication is much faster than division.
368 return GetSizeIJ(level) * (1.0 / kMaxSize);
369}
370
371inline bool S2CellId::is_leaf() const {
372 return int(id_) & 1;
373}
374
375inline bool S2CellId::is_face() const {
376 return (id_ & (lsb_for_level(0) - 1)) == 0;
377}
378
379inline int S2CellId::child_position(int level) const {
380 DCHECK(is_valid());
381 return static_cast<int>(id_ >> (2 * (kMaxLevel - level) + 1)) & 3;
382}
383
384inline S2CellId S2CellId::range_min() const {
385 return S2CellId(id_ - (lsb() - 1));
386}
387
388inline S2CellId S2CellId::range_max() const {
389 return S2CellId(id_ + (lsb() - 1));
390}
391
392inline bool S2CellId::contains(S2CellId const& other) const {
393 DCHECK(is_valid());
394 DCHECK(other.is_valid());
395 return other >= range_min() && other <= range_max();
396}
397
398inline bool S2CellId::intersects(S2CellId const& other) const {
399 DCHECK(is_valid());
400 DCHECK(other.is_valid());
401 return other.range_min() <= range_max() && other.range_max() >= range_min();
402}
403
404inline S2CellId S2CellId::parent(int level) const {
405 DCHECK(is_valid());
406 DCHECK_GE(level, 0);
407 DCHECK_LE(level, this->level());
408 uint64 new_lsb = lsb_for_level(level);
409 return S2CellId((id_ & -new_lsb) | new_lsb);
410}
411
412inline S2CellId S2CellId::parent() const {
413 DCHECK(is_valid());
414 DCHECK(!is_face());
415 uint64 new_lsb = lsb() << 2;
416 return S2CellId((id_ & -new_lsb) | new_lsb);
417}
418
419inline S2CellId S2CellId::child(int position) const {
420 DCHECK(is_valid());
421 DCHECK(!is_leaf());
422 // To change the level, we need to move the least-significant bit two
423 // positions downward. We do this by subtracting (4 * new_lsb) and adding
424 // new_lsb. Then to advance to the given child cell, we add
425 // (2 * position * new_lsb).
426 uint64 new_lsb = lsb() >> 2;
427 return S2CellId(id_ + (2 * position + 1 - 4) * new_lsb);
428}
429
430inline S2CellId S2CellId::child_begin() const {
431 DCHECK(is_valid());
432 DCHECK(!is_leaf());
433 uint64 old_lsb = lsb();
434 return S2CellId(id_ - old_lsb + (old_lsb >> 2));
435}
436
437inline S2CellId S2CellId::child_begin(int level) const {
438 DCHECK(is_valid());
439 DCHECK_GE(level, this->level());
440 DCHECK_LE(level, kMaxLevel);
441 return S2CellId(id_ - lsb() + lsb_for_level(level));
442}
443
444inline S2CellId S2CellId::child_end() const {
445 DCHECK(is_valid());
446 DCHECK(!is_leaf());
447 uint64 old_lsb = lsb();
448 return S2CellId(id_ + old_lsb + (old_lsb >> 2));
449}
450
451inline S2CellId S2CellId::child_end(int level) const {
452 DCHECK(is_valid());
453 DCHECK_GE(level, this->level());
454 DCHECK_LE(level, kMaxLevel);
455 return S2CellId(id_ + lsb() + lsb_for_level(level));
456}
457
458inline S2CellId S2CellId::next() const {
459 return S2CellId(id_ + (lsb() << 1));
460}
461
462inline S2CellId S2CellId::prev() const {
463 return S2CellId(id_ - (lsb() << 1));
464}
465
466inline S2CellId S2CellId::next_wrap() const {
467 DCHECK(is_valid());
468 S2CellId n = next();
469 if (n.id_ < kWrapOffset) return n;
470 return S2CellId(n.id_ - kWrapOffset);
471}
472
473inline S2CellId S2CellId::prev_wrap() const {
474 DCHECK(is_valid());
475 S2CellId p = prev();
476 if (p.id_ < kWrapOffset) return p;
477 return S2CellId(p.id_ + kWrapOffset);
478}
479
480inline S2CellId S2CellId::Begin(int level) {
481 return FromFacePosLevel(0, 0, 0).child_begin(level);
482}
483
484inline S2CellId S2CellId::End(int level) {
485 return FromFacePosLevel(5, 0, 0).child_end(level);
486}
487
488ostream& operator<<(ostream& os, S2CellId const& id);
489
490#ifndef SWIG
491#include<hash_set>
492namespace __gnu_cxx {
493
494
495template<> struct hash<S2CellId> {
496 size_t operator()(S2CellId const& id) const {
497 return static_cast<size_t>(id.id() >> 32) + static_cast<size_t>(id.id());
498 }
499};
500
501
502} // namespace __gnu_cxx
503
504#endif // SWIG
505
506#endif // UTIL_GEOMETRY_S2CELLID_H_
507