diff --git a/docs/.gitignore b/docs/.gitignore new file mode 100644 index 00000000..b439b05c --- /dev/null +++ b/docs/.gitignore @@ -0,0 +1,6 @@ +_site +.sass-cache +.jekyll-metadata +.DS_Store +Gemfile.lock +*~ diff --git a/docs/CNAME b/docs/CNAME new file mode 100644 index 00000000..34fdd53e --- /dev/null +++ b/docs/CNAME @@ -0,0 +1 @@ +s2geometry.io \ No newline at end of file diff --git a/docs/Gemfile b/docs/Gemfile new file mode 100644 index 00000000..021c9a6f --- /dev/null +++ b/docs/Gemfile @@ -0,0 +1,28 @@ +source "https://rubygems.org" +ruby RUBY_VERSION + +# Hello! This is where you manage which Jekyll version is used to run. +# When you want to use a different version, change it below, save the +# file and run `bundle install`. Run Jekyll with `bundle exec`, like so: +# +# bundle exec jekyll serve +# +# This will help ensure the proper Jekyll version is running. +# Happy Jekylling! +# gem "jekyll", "3.4.3" + +# This is the default theme for new Jekyll sites. You may change this to anything you like. +gem "minima" + +# If you want to use GitHub Pages, remove the "gem "jekyll"" above and +# uncomment the line below. To upgrade, run `bundle update github-pages`. +gem "github-pages", group: :jekyll_plugins + +# If you have any plugins, put them here! +group :jekyll_plugins do + gem "jekyll-feed", "~> 0.6" +end + +# Windows does not include zoneinfo files, so bundle the tzinfo-data gem +gem 'tzinfo-data', platforms: [:mingw, :mswin, :x64_mingw, :jruby] + diff --git a/docs/LICENSE b/docs/LICENSE new file mode 100644 index 00000000..b014ed5b --- /dev/null +++ b/docs/LICENSE @@ -0,0 +1,204 @@ + Apache License + Version 2.0, January 2004 + http://www.apache.org/licenses/ + + TERMS AND CONDITIONS FOR USE, REPRODUCTION, AND DISTRIBUTION + + 1. Definitions. + + "License" shall mean the terms and conditions for use, reproduction, + and distribution as defined by Sections 1 through 9 of this document. + + "Licensor" shall mean the copyright owner or entity authorized by + the copyright owner that is granting the License. + + "Legal Entity" shall mean the union of the acting entity and all + other entities that control, are controlled by, or are under common + control with that entity. For the purposes of this definition, + "control" means (i) the power, direct or indirect, to cause the + direction or management of such entity, whether by contract or + otherwise, or (ii) ownership of fifty percent (50%) or more of the + outstanding shares, or (iii) beneficial ownership of such entity. + + "You" (or "Your") shall mean an individual or Legal Entity + exercising permissions granted by this License. + + "Source" form shall mean the preferred form for making modifications, + including but not limited to software source code, documentation + source, and configuration files. + + "Object" form shall mean any form resulting from mechanical + transformation or translation of a Source form, including but + not limited to compiled object code, generated documentation, + and conversions to other media types. + + "Work" shall mean the work of authorship, whether in Source or + Object form, made available under the License, as indicated by a + copyright notice that is included in or attached to the work + (an example is provided in the Appendix below). + + "Derivative Works" shall mean any work, whether in Source or Object + form, that is based on (or derived from) the Work and for which the + editorial revisions, annotations, elaborations, or other modifications + represent, as a whole, an original work of authorship. For the purposes + of this License, Derivative Works shall not include works that remain + separable from, or merely link (or bind by name) to the interfaces of, + the Work and Derivative Works thereof. + + "Contribution" shall mean any work of authorship, including + the original version of the Work and any modifications or additions + to that Work or Derivative Works thereof, that is intentionally + submitted to Licensor for inclusion in the Work by the copyright owner + or by an individual or Legal Entity authorized to submit on behalf of + the copyright owner. For the purposes of this definition, "submitted" + means any form of electronic, verbal, or written communication sent + to the Licensor or its representatives, including but not limited to + communication on electronic mailing lists, source code control systems, + and issue tracking systems that are managed by, or on behalf of, the + Licensor for the purpose of discussing and improving the Work, but + excluding communication that is conspicuously marked or otherwise + designated in writing by the copyright owner as "Not a Contribution." + + "Contributor" shall mean Licensor and any individual or Legal Entity + on behalf of whom a Contribution has been received by Licensor and + subsequently incorporated within the Work. + + 2. Grant of Copyright License. Subject to the terms and conditions of + this License, each Contributor hereby grants to You a perpetual, + worldwide, non-exclusive, no-charge, royalty-free, irrevocable + copyright license to reproduce, prepare Derivative Works of, + publicly display, publicly perform, sublicense, and distribute the + Work and such Derivative Works in Source or Object form. + + 3. Grant of Patent License. Subject to the terms and conditions of + this License, each Contributor hereby grants to You a perpetual, + worldwide, non-exclusive, no-charge, royalty-free, irrevocable + (except as stated in this section) patent license to make, have made, + use, offer to sell, sell, import, and otherwise transfer the Work, + where such license applies only to those patent claims licensable + by such Contributor that are necessarily infringed by their + Contribution(s) alone or by combination of their Contribution(s) + with the Work to which such Contribution(s) was submitted. If You + institute patent litigation against any entity (including a + cross-claim or counterclaim in a lawsuit) alleging that the Work + or a Contribution incorporated within the Work constitutes direct + or contributory patent infringement, then any patent licenses + granted to You under this License for that Work shall terminate + as of the date such litigation is filed. + + 4. Redistribution. You may reproduce and distribute copies of the + Work or Derivative Works thereof in any medium, with or without + modifications, and in Source or Object form, provided that You + meet the following conditions: + + (a) You must give any other recipients of the Work or + Derivative Works a copy of this License; and + + (b) You must cause any modified files to carry prominent notices + stating that You changed the files; and + + (c) You must retain, in the Source form of any Derivative Works + that You distribute, all copyright, patent, trademark, and + attribution notices from the Source form of the Work, + excluding those notices that do not pertain to any part of + the Derivative Works; and + + (d) If the Work includes a "NOTICE" text file as part of its + distribution, then any Derivative Works that You distribute must + include a readable copy of the attribution notices contained + within such NOTICE file, excluding those notices that do not + pertain to any part of the Derivative Works, in at least one + of the following places: within a NOTICE text file distributed + as part of the Derivative Works; within the Source form or + documentation, if provided along with the Derivative Works; or, + within a display generated by the Derivative Works, if and + wherever such third-party notices normally appear. The contents + of the NOTICE file are for informational purposes only and + do not modify the License. You may add Your own attribution + notices within Derivative Works that You distribute, alongside + or as an addendum to the NOTICE text from the Work, provided + that such additional attribution notices cannot be construed + as modifying the License. + + You may add Your own copyright statement to Your modifications and + may provide additional or different license terms and conditions + for use, reproduction, or distribution of Your modifications, or + for any such Derivative Works as a whole, provided Your use, + reproduction, and distribution of the Work otherwise complies with + the conditions stated in this License. + + 5. Submission of Contributions. Unless You explicitly state otherwise, + any Contribution intentionally submitted for inclusion in the Work + by You to the Licensor shall be under the terms and conditions of + this License, without any additional terms or conditions. + Notwithstanding the above, nothing herein shall supersede or modify + the terms of any separate license agreement you may have executed + with Licensor regarding such Contributions. + + 6. Trademarks. This License does not grant permission to use the trade + names, trademarks, service marks, or product names of the Licensor, + except as required for reasonable and customary use in describing the + origin of the Work and reproducing the content of the NOTICE file. + + 7. Disclaimer of Warranty. Unless required by applicable law or + agreed to in writing, Licensor provides the Work (and each + Contributor provides its Contributions) on an "AS IS" BASIS, + WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or + implied, including, without limitation, any warranties or conditions + of TITLE, NON-INFRINGEMENT, MERCHANTABILITY, or FITNESS FOR A + PARTICULAR PURPOSE. You are solely responsible for determining the + appropriateness of using or redistributing the Work and assume any + risks associated with Your exercise of permissions under this License. + + 8. Limitation of Liability. In no event and under no legal theory, + whether in tort (including negligence), contract, or otherwise, + unless required by applicable law (such as deliberate and grossly + negligent acts) or agreed to in writing, shall any Contributor be + liable to You for damages, including any direct, indirect, special, + incidental, or consequential damages of any character arising as a + result of this License or out of the use or inability to use the + Work (including but not limited to damages for loss of goodwill, + work stoppage, computer failure or malfunction, or any and all + other commercial damages or losses), even if such Contributor + has been advised of the possibility of such damages. + + 9. Accepting Warranty or Additional Liability. While redistributing + the Work or Derivative Works thereof, You may choose to offer, + and charge a fee for, acceptance of support, warranty, indemnity, + or other liability obligations and/or rights consistent with this + License. However, in accepting such obligations, You may act only + on Your own behalf and on Your sole responsibility, not on behalf + of any other Contributor, and only if You agree to indemnify, + defend, and hold each Contributor harmless for any liability + incurred by, or claims asserted against, such Contributor by reason + of your accepting any such warranty or additional liability. + + END OF TERMS AND CONDITIONS + + APPENDIX: How to apply the Apache License to your work. + + To apply the Apache License to your work, attach the following + boilerplate notice, with the fields enclosed by brackets "[]" + replaced with your own identifying information. (Don't include + the brackets!) The text should be enclosed in the appropriate + comment syntax for the file format. We also recommend that a + file or class name and description of purpose be included on the + same "printed page" as the copyright notice for easier + identification within third-party archives. + + Copyright [yyyy] [name of copyright owner] + + Licensed under the Apache License, Version 2.0 (the "License"); + you may not use this file except in compliance with the License. + You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + + Unless required by applicable law or agreed to in writing, software + distributed under the License is distributed on an "AS IS" BASIS, + WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + See the License for the specific language governing permissions and + limitations under the License. + + +© 2017 GitHub, Inc. \ No newline at end of file diff --git a/docs/README.md b/docs/README.md new file mode 100644 index 00000000..5d00132b --- /dev/null +++ b/docs/README.md @@ -0,0 +1,4 @@ +# s2geometry.github.io + +Website for S2 Geometry project + diff --git a/docs/_config.yml b/docs/_config.yml new file mode 100644 index 00000000..e8787ffe --- /dev/null +++ b/docs/_config.yml @@ -0,0 +1,57 @@ +# Welcome to Jekyll! +# +# This config file is meant for settings that affect your whole blog, values +# which you are expected to set up once and rarely edit after that. If you find +# yourself editing this file very often, consider using Jekyll's data files +# feature for the data you need to update frequently. +# +# For technical reasons, this file is *NOT* reloaded automatically when you use +# 'bundle exec jekyll serve'. If you change this file, please restart the server process. + +# Site settings +# These are used to personalize your new site. If you look in the HTML files, +# you will see them accessed via {{ site.title }}, {{ site.email }}, and so on. +# You can create any custom variable you would like, and they will be accessible +# in the templates via {{ site.myvariable }}. + +# Not currently used +title: S2Geometry +email: s2geometry-io@googlegroups.com +description: The s2geometry.io website + +# The base hostname & protocol for your site, e.g. http://example.com +# This site.url path is required by the RSS feed.xml +url: "http://s2geometry.io" + +# the subpath of your site. The main Abseil website is not on a subpath +baseurl: "" + +# host: 0.0.0.0 is required for local development. +host: 0.0.0.0 + +# Build settings +markdown: kramdown +theme: minima + +# Minima Navigation Config +header_pages: + - getS2.md + - about/index.md + - devguide/index.md + - community/index.md + - resources/index.md + +# RSS / Blog configuration settings +gems: + - jekyll-feed +timezone: "America/New_York" + +# Set these false for production, true for staging server +future: false +unpublished: false +show_drafts: false + +# Exclusions +exclude: + - Gemfile + - Gemfile.lock diff --git a/docs/_includes/header.html b/docs/_includes/header.html new file mode 100644 index 00000000..72e2d9d8 --- /dev/null +++ b/docs/_includes/header.html @@ -0,0 +1,18 @@ + diff --git a/docs/about/index.md b/docs/about/index.md new file mode 100644 index 00000000..e7d83fa8 --- /dev/null +++ b/docs/about/index.md @@ -0,0 +1,29 @@ +--- +title: About S2 +--- + +This guide provides background information for starting work with +the S2 Geometry Library. + +* The S2 [Platforms Guide](platforms) denotes the prerequisites + for installation of S2 on your system. This guide also + provides build instructions for S2. +* The S2 [Overview](overview) provides a conceptual introduction + to S2, and should be read before diving into the developer + guide. + +S2 consists of source code repositories for C++ and Java. Upon +completing the overview, you will probably want to run through +the Quickstart: + +* [C++ Quickstart](/devguide/cpp/quickstart) + +(A Java Quickstart is envisioned soon.) + +Ready for More? + +* Read about [S2 Cells](/devguide/s2cell_hierarchy) +* Learn about the [Basic Types](/devguide/basic_types) +* Understand the main [S2ShapeIndex](/devguide/s2shapeindex) + Class +* Discover [Finding Nearby Edges](/devguide/s2closestedgequery) diff --git a/docs/about/install.md b/docs/about/install.md new file mode 100644 index 00000000..3185b99a --- /dev/null +++ b/docs/about/install.md @@ -0,0 +1,247 @@ +--- +title: S2 Installation +--- + +This document lists the platform requirements and installation +instructions for working with the S2 Geometry Library in both +C++ and Python. (The S2 Python interfaces uses SWIG to interact +with the C++ code.) + +## Requirements + +Running the S2 code requires the following: + +* A MacOSX or Linux platform. (Windows is not supported at this time.) +* POSIX support (for `getrusage()`). +* A compatible C++ compiler *supporting at least C++11*, such as + [g++](https://gcc.gnu.org/){:target="_blank"} >= 4.7. +* Installation of the following libraries: + * [OpenSSL](https://github.com/openssl/openssl){:target="_blank"} (for its + bignum library) + * [gflags command line flags](https://github.com/gflags/gflags) + {:target="_blank"} (optional, disabled by default) + * [glog logging module](https://github.com/google/glog) {:target="_blank"} + (optional, disabled by default) + * [googletest testing framework](https://github.com/google/googletest) + {:target="_blank"} + * [SWIG](https://github.com/swig/swig) (optional, for Python interface) +* [Git](https://git-scm.com/) for interacting with the S2 source code + repository, which is contained on [GitHub](http://github.com). To install + Git, consult the [Set Up Git](https://help.github.com/articles/set-up-git/) + guide on GitHub. + +Installation instructions for the above components are [noted below](#Install). + +

+Note: this installation guide uses CMake as the official build system for +S2, which is supported on most major platforms and compilers. The S2 +source code assumes you are using CMake and contains +CMakeList.txt files for that purpose. +

+ +Although you are free to use your own build system, most of the documentation +within this guide will assume you are using +[CMake](https://cmake.org/){:target="_blank"}. + +## Installation Instructions {#Install} + +Note: thorough testing has only been done on Ubuntu 14.04.3 +and MacOSX 10.12. + +### Setting Up Your Development Environment (Linux) + +Note: we recommend you use a Linux package manager such as +`apt-get`. Alternatively, you may build the dependent +libraries from source. + +1\. Install the additional libraries S2 code requires: + +
+$ sudo apt-get install libgflags-dev libgoogle-glog-dev libgtest-dev libssl-dev
+Reading package lists... Done
+Building dependency tree
+Reading state information... Done
+...
+After this operation, 3,090 kB of additional disk space will be used.
+Do you want to continue? [Y/n] Y
+...
+Setting up libgtest-dev (1.6.0-1ubuntu6) ...
+Processing triggers for libc-bin (2.19-0ubuntu6.13) ...
+$
+
+ +2\. Install CMake + +
+$ sudo apt-get install cmake
+Reading package lists... Done
+Building dependency tree
+Reading state information... Done
+...
+After this operation, 16.6 MB of additional disk space will be used.
+Do you want to continue? [Y/n] Y
+...
+Setting up cmake (2.8.12.2-0ubuntu3) ...
+$
+
+ +### Setting Up Your Development Environment (MacOSX) + +Note: we recommend you use a MacOSX package manager such as +[MacPorts](https://macports.org) or [Homebrew](https://brew.sh). +Alternatively, you may build the dependent libraries from source. + +1\. Install the additional libraries S2 code requires: + +
+# MacPorts
+$ sudo port install gflags google-glog openssl
+
+# Homebrew
+$ sudo brew install gflags glog openssl
+
+ +2\. Download Googletest +[release 1.8.0](https://github.com/google/googletest/releases/tag/release-1.8.0){:target="_blank"} +and unpack it in a directory of your choosing. (Take note of this +directory as you will need to point CMake to it.) + +3\. Install CMake + +
+# Note: XCode requires command-line tools, which can be installed with:
+$ xcode-select --install
+
+# MacPorts
+$ sudo port install cmake
+
+# Homebrew
+$ sudo brew install cmake
+
+ +### Getting the S2 Code + +The [Reference implementation of S2](https://github.com/google/s2geometry) +is written in C++. Once you have CMake and Git installed, you can obtain +the S2 code from its repository on GitHub: + +
+# Change to the directory where you want to create the code repository
+$ cd ~
+$ mkdir Source; cd Source
+
+ +Clone the S2 code into your development directory: + +
+$ git clone https://github.com/google/s2geometry.git
+Cloning into 's2geometry'...
+remote: Counting objects: 5672, done.
+remote: Compressing objects: 100% (52/52), done.
+remote: Total 5672 (delta 21), reused 36 (delta 17), pack-reused 5601
+Receiving objects: 100% (5672/5672), 7.15 MiB | 27.21 MiB/s, done.
+Resolving deltas: 100% (4503/4503), done.
+$
+
+ +Git will create the repository within a directory named `s2geometry`. +Navigate into this directory. + +

+Alternatively, you can download an archive of the S2 library as a +ZIP file +and unzip it within your development directory. + +

+$ unzip path_to_ZIP_file/s2geometry-master.zip
+
+

+ +### Compiling S2 + +Once you have installed S2's requirements and the S2 library, you +are ready to build S2. + +1\. Configure the make files using CMake: + +
+$ cd s2geometry  # or geometry-master if you installed from the ZIP file
+$ mkdir build
+$ cd build
+# See Notes below.
+$ cmake -DWITH_GFLAGS=ON -WITH_GTEST=ON \
+-DGTEST_ROOT=googletest_root_dir \
+-DOPENSSL_INCLUDE_DIR=openssl_include_dir ..
+-- Found OpenSSL: /usr/lib/libcrypto.dylib (found version "1.0.2m")
+-- Looking for pthread.h
+-- Looking for pthread.h - found
+-- Looking for pthread_create
+-- Looking for pthread_create - found
+-- Found Threads: TRUE
+-- Could NOT find SWIG (missing: SWIG_EXECUTABLE SWIG_DIR)
+-- Found PythonInterp: /Library/Frameworks/Python.framework/Versions/2.7/bin/python (found version "2.7.11")
+-- Found PythonLibs: /Library/Frameworks/Python.framework/Versions/2.7/lib/libpython2.7.dylib (found version "2.7.11")
+GTEST_ROOT: /Users/shreck/SOURCE/googletest
+-- Configuring done
+-- Generating done
+-- Build files have been written to: /Users/shreck/SOURCE/s2geometry-master/build
+$
+
+ +Notes: + +* `-DWITH_GFLAGS` and `-DWITH_GTEST` may be omitted to avoid depending + on those packages. +* `-DGTEST_ROOT` and `-DOPENSSL_INCLUDE_DIR` must be absolute paths. +* `-DGTEST_ROOT` is usually `/usr/src/gtest` on Linux systems; on + MacOSX use the directory where you installed GoogleTest. +* MacOSX/Homebrew users may need to set `-DOPENSSL_INCLUDE_DIR` to + enable CMake to find your openssl `include` directory. +* You can omit the `DGTEST_ROOT` flag to skip tests. + +2\. Make the S2 binary (and associated tests if `GTEST_ROOT` was + specified in CMake: + +
+$ make
+Scanning dependencies of target gtest
+...
+Scanning dependencies of target s2
+...
+[ 35%] Built target s2
+...
+[100%] Built target point_index
+$
+
+ +3\. Run the S2 tests (if `GTEST_ROOT` was specified in CMake): + +
+$ make test
+Running tests...
+Test project /.../s2geometry-master/build
+      Start  1: id_set_lexicon_test
+...
+83/83 Test #83: value_lexicon_test .............................   Passed    0.01 sec
+100% tests passed, 0 tests failed out of 83
+Total Test time (real) =  78.67 sec
+$
+
+ + +4\. Install the built S2 library: + +
+$ make install
+[  1%] Built target gtest
+...
+[100%] Built target point_index
+Install the project...
+-- Install configuration: ""
+...
+$
+
+ +Congratulations! You've successfully installed S2, built the library, and run +all of the tests! You're now ready to move on to the +[C++ Quickstart](/devguide/cpp/quickstart). diff --git a/docs/about/overview.md b/docs/about/overview.md new file mode 100644 index 00000000..b9368ccc --- /dev/null +++ b/docs/about/overview.md @@ -0,0 +1,317 @@ +--- +title: Overview +--- + +S2 is a library for spherical geometry that aims to have the same +robustness, flexibility, and performance as the very best planar +geometry libraries. + +## Spherical Geometry + +Let's break down the elements of this goal. First, why spherical +geometry? (By the way, the name "S2" is derived from the mathematical +notation for the unit sphere, *S²*.) + +Traditional cartography is based on *map projections*, which are simply +functions that map points on the Earth's surface to points on a planar +map. Map projections create distortions due to the fact that the shape +of the Earth is not very close to the shape of a plane. For example, +the well-known Mercator projection is discontinuous along the 180 degree +meridian, has large scale distortions at high latitudes, and cannot +represent the north and south poles at all. Other projections make +different compromises, but no planar projection does a good job of +representing the entire surface of the Earth. + +S2 approaches this problem by working exclusively with *spherical +projections*. As the name implies, spherical projections map points on +the Earth's surface to a perfect mathematical sphere. Such mappings +still create some distortion, of course, because the Earth is not quite +spherical--but as it turns out, the Earth is much closer to being a +sphere than a plane. With spherical projections, it is possible to +approximate the entire Earth's surface with a maximum distortion of +0.56%. Perhaps more importantly, spherical projections preserve the +correct topology of the Earth -- there are no singularities or +discontinuities to deal with. + +Why not project onto an ellipsoid? (The Earth isn't quite ellipsoidal +either, but it is even closer to being an ellipsoid than a sphere.) +The answer relates to the other goals stated above, namely performance +and robustness. Ellipsoidal operations are still orders of magnitude +slower than the corresponding operations on a sphere. Furthermore, +robust geometric algorithms require the implementation of exact +geometric predicates that are not subject to numerical errors. While +this is fairly straightforward for planar geometry, and somewhat harder +for spherical geometry, it is not known how to implement all of the +necessary predicates for ellipsoidal geometry. + +## Robustness + +Which brings up the question, what do we mean by "robust"? + +In the S2 library, the core operations are designed to be 100% robust. +This means that each operation makes strict mathematical guarantees +about its output, and is implemented in such a way that it meets those +guarantees for all possible valid inputs. For example, if you compute +the intersection of two polygons, not only is the output guaranteed to +be topologically correct (up to the creation of degeneracies), but it +is also guaranteed that the boundary of the output stays within a +user-specified tolerance of true, mathematically exact result. + +Robustness is very important when building higher-level algorithms, +since unexpected results from low-level operations can be very +difficult to handle. S2 achieves this goal using a combination of +techniques from computational geometry, including *conservative error +bounds*, *exact geometric predicates*, and *snap rounding*. + +S2 attempts to be precise both in terms of mathematical definitions +(e.g., whether regions include their boundaries, and how degeneracies +are handled) and numerical accuracy (e.g., minimizing cancellation +error). + +## Flexibility + +S2 is organized as a toolkit that gives clients as much control as +possible. For example, S2 makes it easy to implement your own geometry +subtypes that store data in a format of your choice (the `S2Shape` +interface). Similarly, most operations are designed to give clients +fine control over the semantics, e.g. whether polygon boundaries are +considered to be open or closed (or semi-open). S2 has APIs at several +different levels to accommodate various uses of the library. + +## Performance + +S2 is designed to have good performance on large geographic datasets. +Most operations are accelerated using an in-memory edge index data +structure (`S2ShapeIndex`). For example if you have a million +polygons, finding the polygon(s) that contain a given point typically +takes a few hundred nanoseconds. Similarly it is fast to find objects +that are near each other, such as finding all the places of business +near a given road, or all the roads near a given location. + +Many operations on indexed data also have output-sensitive running +times. For example, if you intersect the Pacific Ocean with a small +rectangle near Seattle, the running time essentially depends on the +complexity of the Pacific Ocean near Seattle, not the complexity of the +Pacific Ocean overall. + +## Scope + +The S2 library provides the following: + +* Representations of angles, intervals, latitude-longitude points, + unit vectors, and so on, and various operations on these types. + +* Geometric shapes over the unit sphere, such as spherical caps + ("discs"), latitude-longitude rectangles, polylines, and polygons. + +* Robust constructive operations (e.g., union) and boolean predicates + (e.g., containment) for arbitrary collections of points, polylines, + and polygons. + +* Fast in-memory indexing of collections of points, polylines, and + polygons. + +* Algorithms for measuring distances and finding nearby objects. + +* Robust algorithms for snapping and simplifying geometry (with + accuracy and topology guarantees). + +* A collection of efficient yet exact mathematical predicates for + testing relationships among geometric objects. + +* Support for spatial indexing, including the ability to approximate regions + as collections of discrete "S2 cells". This feature makes it easy to + build large distributed spatial indexes. + +* Extensive testing on Google's vast collection of geographic data. + +* Flexible Apache 2.0 license. + +On the other hand, the following are outside the scope of S2: + +* Planar geometry. (There are many fine existing planar geometry + libraries to choose from.) + +* Conversions to/from common GIS formats. (To read such formats, use + an external library such as + [OGR](http://gdal.org/1.11/ogr/){:target="_blank"}.) + +## Language Support + +The [reference implementation of the S2 +library](https://github.com/google/s2geometry) is written in C++. Versions in +other languages are ported from the C++ code, and may not have the same +robustness, performance, or features as the C++ version. As of November 2017, +the "official" ports are as follows: + +* [S2 library in Go](https://github.com/golang/geo). This project is being + actively developed and aims to have full equivalence with the C++ version. + +* [S2 library in Java](https://github.com/google/s2-geometry-library-java). + Be aware that this port is based on the 2011 release of the S2 library, and + does not have the same robustness, performance, or features as the current + C++ implementation. + +* [S2 library in Python](https://github.com/google/s2geometry/tree/master/src/python). + A subset of the library has been ported to Python using SWIG, and is + included in the `python` directory of the C++ repository. The supported + feature set includes polylines, polygons, discs, rectangles, and `S2Cell` + coverings. + +## Layered Design + +The S2 library is structured as a toolkit with various layers. At the +lowest level, it provides a set of operations on points, edges, and simple +shapes such as rectangles and discs. For example, this includes classes +such as `S2Point`, `S2CellId`, `S1Angle`, `S1Interval`, `S2LatLng`, +`S2Region`. It also includes the exact predicates defined in +`s2predicates.h` and `s2edge_crossings.h`, utility functions such as those +in `s2edge_distances.h`, and classes related to the `S2Cell` hierarchy (a +hierarchical decomposition of the sphere) such as `S2CellId`, `S2Region`, +and `S2RegionCoverer`. + +The next layer provides flexible support for polygonal geometry, including +point sets and polylines. It is built around the following core classes: + +* `S2Builder`: a robust tool for constructing, and simplifying + geometry, built on the framework of snap rounding. + +* `S2ShapeIndex`: an indexed collection of points, polylines, and/or + polygons, possibly overlapping. + +* `S2BooleanOperation`: evaluates boolean operations and predicates for + arbitrary collections of polygonal geometry. + +* `S2ClosestEdgeQuery`: measures distances and finds closest edges for + arbitrary collections of polygonal geometry. + +These are only the most important classes, but what they have in common +is that (1) they work with arbitrary collections of points, polylines, +and polygons; (2) they allow clients to control the underlying geometry +representation (via the `S2Shape` interface), and (3) they give clients +fine control over the semantics of operations (e.g., how to handle +degeneracies, whether polygons are open or closed, whether and how the +output should be snapped, etc). + +The final layer consists of convenience classes, such as `S2Polygon`, +`S2Loop`, and `S2Polyline`. These types have wide interfaces and +convenience methods that implement many of the operations above (e.g., +intersection, distances), however they don't offer the same flexibility +in terms of data representation and control over semantics (e.g., +polygon boundaries are always semi-open). + +For example, you can measure the distance from a point to a polygon as +follows: + + S2Point p = ...; + S2Polygon a = ...; + S1Angle dist = a.GetDistance(p); + +whereas with the `S2ShapeIndex` layer you would write + + S2Point p = ...; + S2Polygon a = ...; + S2ShapeIndex index; + index.Add(absl::make_unique(&a)); + S2ClosestEdgeQuery query(&index); + S1ChordAngle dist = query.GetDistance(p); + +The first case is simpler, but the second case supports many more +options (e.g. measuring the distance to a collection of geometry, or +finding the objects within a certain distance). + +## Design Choices + +This section summarizes some of the major design choices made by the library. + +### Spherical Geodesic Edges + +In S2, all edges are "spherical geodesics", i.e. shortest paths on the +sphere. For example, the edge between two points 100 meters apart on +opposite side of the North pole goes directly through the North pole. +In contrast, the same edge in the standard "plate carrée" +projection (i.e., raw latitude/longitude coordinates) would follow a +line of constant latitude, which in this case happens to be a +semicircular path that always stays exactly 50 meters away from the +North pole. (Lines of latitude except for the equator are never +shortest paths.) + +With geodesic edges there are no special cases near the poles or the +180 degree meridian; for example, if two points 10km apart are +separated by the 180 degree meridian, the edge between them simply +crosses the meridian rather than going all the way around the other +side of the Earth. + +### Orientation Matters + +In S2, polygon loops contain the region on their left. This differs +from planar libraries that might define loops as being "clockwise" or +"counter-clockwise", or that accept loops of either orientation. None +of these concepts really applies to polygons drawn on a sphere; for +example, consider a loop that goes around the equator. Is it clockwise +or counter-clockwise? (For those familiar with the term, the concept +of "winding number" does not exist on the sphere.) + +For these reasons, S2 follows the "interior is on the left" rule. (This +is the same as the counter-clockwise rule for small loops, but differs +for very large loops that enclose more than half of the sphere.) + +### Unit Vectors vs. Latitude/Longitude Pairs + +In S2, points are represented internally as unit-length vectors (points +on the surface of a three-dimensional unit sphere) as opposed to +traditional (latitude, longitude) pairs. This is for two reasons: + +* Unit vectors are much more efficient when working with geodesic + edges. Using (latitude, longitude) pairs would require constantly + evaluating trigonometric functions (sin, cos, etc), which is slow + even on modern CPU architectures. + +* Unit vectors allow exact geometric predicates to be evaluated + efficiently. To do this with (latitude, longitude) pairs, you + would essentially need to convert them to the unit vector + representation first. + +For precision and efficiency, coordinates are representated as +double-precision numbers. Robustness is achieved through the use of +the techniques mentioned previously: conservative error bounds (which +measure the error in certain calculations), exact geometric predicates +(which determine the true mathematical relationship among geometric +objects), and snap rounding (a technique that allows rounding results +to finite precision without creating topological errors). + +### Classes vs. Functions + +Most algorithms in S2 are implemented as classes rather than functions. +So for example, rather than having a `GetDistance` function, there is +an `S2ClosestEdgeQuery` class with a `GetDistance` method. This design +has two advantages: + +* It allows caching and reuse of data structures. For example, if an + algorithm needs to make thousands of distance queries, it can + declare one `S2ClosestEdgeQuery` object and call its methods + thousands of times. This provides the opportunity to avoid + recomputing data that can be used from one query to the next. + +* It provides a convenient way to specify options (via a nested + `Options` class), and to avoid respecifying those options many + times when repeated operations need to be performed. + +### Geocentric vs. Geodetic Coordinates + +Geocentric and geodetic coordinates are two methods for defining a +(latitude, longitude) coordinate system over the Earth's surface. The +question sometimes arises: which system does the S2 library use? + +The answer is: neither (or both). Points in S2 are modeled as lying on +a perfect sphere (just as points in a traditional planar geometry +library are modeled as lying on a perfect plane). How points are +mapped onto the sphere is outside the scope of S2; it works equally +well with geocentric or geodetic coordinates. (Actual geographic +datasets use geodetic coordinates exclusively.) + +## Authors + +The S2 library was written primarily by Eric Veach. Other significant +contributors include Jesse Rosenstock, Eric Engle (Java port lead), Robert +Snedegar (Go port lead), Julien Basch, and Tom Manshreck diff --git a/docs/about/platforms.md b/docs/about/platforms.md new file mode 100644 index 00000000..33393194 --- /dev/null +++ b/docs/about/platforms.md @@ -0,0 +1,273 @@ +--- +title: S2 Installation +--- + +This document lists the platform requirements and installation +instructions for working with the S2 Geometry Library in both +C++ and Python. (The S2 Python interfaces uses SWIG to interact +with the C++ code.) + +## Requirements + +Running the S2 code requires the following: + +* A MacOSX or Linux platform. (Windows is not supported at this time.) +* POSIX support (for `getrusage()`). +* A compatible C++ compiler *supporting at least C++11*, such as + [g++](https://gcc.gnu.org/){:target="_blank"} >= 4.7. +* Installation of the following libraries: + * [OpenSSL](https://github.com/openssl/openssl){:target="_blank"} (for its + bignum library) + * [gflags command line flags](https://github.com/gflags/gflags) + {:target="_blank"} (optional, disabled by default) + * [glog logging module](https://github.com/google/glog) {:target="_blank"} + (optional, disabled by default) + * [googletest testing framework](https://github.com/google/googletest) + {:target="_blank"} + * [SWIG](https://github.com/swig/swig) (optional, for Python interface) +* [Git](https://git-scm.com/) for interacting with the S2 source code + repository, which is contained on [GitHub](http://github.com). To install + Git, consult the [Set Up Git](https://help.github.com/articles/set-up-git/) + guide on GitHub. + +Installation instructions for the above components are [noted below](#install). + +

+Note: this installation guide uses CMake as the official build system for +S2, which is supported on most major platforms and compilers. The S2 +source code assumes you are using CMake and contains +CMakeList.txt files for that purpose. +

+ +Although you are free to use your own build system, most of the documentation +within this guide will assume you are using +[CMake](https://cmake.org/){:target="_blank"}. + +## Installation Instructions {#install} + +Note: thorough testing has only been done on Ubuntu 14.04.3 +and MacOSX 10.12. + +### Setting Up Your Development Environment (Linux) + +Note: we recommend you use a Linux package manager such as +`apt-get`. Alternatively, you may build the dependent +libraries from source. + +1\. Install the additional libraries S2 code requires: + +
+$ sudo apt-get install libgflags-dev libgoogle-glog-dev libgtest-dev libssl-dev
+Reading package lists... Done
+Building dependency tree
+Reading state information... Done
+...
+After this operation, 3,090 kB of additional disk space will be used.
+Do you want to continue? [Y/n] Y
+...
+Setting up libgtest-dev (1.6.0-1ubuntu6) ...
+Processing triggers for libc-bin (2.19-0ubuntu6.13) ...
+$
+
+ +2\. (Optional) Install SWIG + +If you will be developing in Python, you should also install the +[SWIG](http://www.swig.org){:target="_blank"} library). Most +Python code will invoke the C++ code through this library. + +
+$ sudo apt-get install swig
+$
+
+ +3\. Install CMake + +
+$ sudo apt-get install cmake
+Reading package lists... Done
+Building dependency tree
+Reading state information... Done
+...
+After this operation, 16.6 MB of additional disk space will be used.
+Do you want to continue? [Y/n] Y
+...
+Setting up cmake (2.8.12.2-0ubuntu3) ...
+$
+
+ +### Setting Up Your Development Environment (MacOSX) + +Note: we recommend you use a MacOSX package manager such as +[MacPorts](https://macports.org) or [Homebrew](https://brew.sh). +Alternatively, you may build the dependent libraries from source. + +1\. Install the additional libraries S2 code requires: + +
+# MacPorts
+$ sudo port install gflags google-glog openssl
+
+# Homebrew
+$ brew install gflags glog openssl
+
+ +2\. Download Googletest +[release 1.8.0](https://github.com/google/googletest/releases/tag/release-1.8.0){:target="_blank"} +and unpack it in a directory of your choosing. (Take note of this +directory as you will need to point CMake to it.) + +3\. (Optional) Install SWIG + +If you will be developing in Python, you should also install the +[SWIG](http://www.swig.org){:target="_blank"} library). Most Python code will +invoke the C++ code through this library. + +
+# MacPorts
+$ sudo port install swig
+
+# Homebrew
+$ brew install swig
+
+ +4\. Install CMake + +
+# Note: XCode requires command-line tools, which can be installed with:
+$ xcode-select --install
+
+# MacPorts
+$ sudo port install cmake
+
+# Homebrew
+$ brew install cmake
+
+ +### Getting the S2 Code + +The [Reference implementation of S2](https://github.com/google/s2geometry) +is written in C++. Once you have CMake and Git installed, you can obtain +the S2 code from its repository on GitHub: + +
+# Change to the directory where you want to create the code repository
+$ cd ~
+$ mkdir Source; cd Source
+
+ +Clone the S2 code into your development directory: + +
+$ git clone https://github.com/google/s2geometry.git
+Cloning into 's2geometry'...
+remote: Counting objects: 5672, done.
+remote: Compressing objects: 100% (52/52), done.
+remote: Total 5672 (delta 21), reused 36 (delta 17), pack-reused 5601
+Receiving objects: 100% (5672/5672), 7.15 MiB | 27.21 MiB/s, done.
+Resolving deltas: 100% (4503/4503), done.
+$
+
+ +Git will create the repository within a directory named `s2geometry`. +Navigate into this directory. + +

+Alternatively, you can download an archive of the S2 library as a +ZIP file +and unzip it within your development directory. + +

+$ unzip path_to_ZIP_file/s2geometry-master.zip
+
+

+ +### Compiling S2 + +Once you have installed S2's requirements and the S2 library, you +are ready to build S2. + +1\. Configure the make files using CMake: + +
+$ cd s2geometry  # or geometry-master if you installed from the ZIP file
+$ mkdir build
+$ cd build
+# See Notes below on whether/when to use these cmake args
+$ cmake -DWITH_GFLAGS=ON -WITH_GTEST=ON \
+-DGTEST_ROOT=googletest_root_dir \
+-DOPENSSL_INCLUDE_DIR=openssl_include_dir ..
+-- Found OpenSSL: /usr/lib/libcrypto.dylib (found version "1.0.2m")
+-- Looking for pthread.h
+-- Looking for pthread.h - found
+-- Looking for pthread_create
+-- Looking for pthread_create - found
+-- Found Threads: TRUE
+-- Found PythonInterp: /Library/Frameworks/Python.framework/Versions/2.7/bin/python (found version "2.7.11")
+-- Found PythonLibs: /Library/Frameworks/Python.framework/Versions/2.7/lib/libpython2.7.dylib (found version "2.7.11")
+GTEST_ROOT: /absolute_path/googletest
+-- Configuring done
+-- Generating done
+-- Build files have been written to: /absolute_path/s2geometry-master/build
+$
+
+ +Notes: + +* `-DWITH_GFLAGS` and `-DWITH_GTEST` may be omitted to avoid depending + on those packages. +* `-DGTEST_ROOT` and `-DOPENSSL_INCLUDE_DIR` must be absolute paths. +* `-DGTEST_ROOT`is the root directory for GoogleTest (e.g. `/usr/src/gtest` on + Linux systems, or the directory you directly installed GoogleTest within, + such as /Users/username/source/googletest on MacOSX). +* MacOSX/Homebrew users may need to set `-DOPENSSL_INCLUDE_DIR` to enable + CMake to find your openssl `include` directory (e.g. + `/usr/local/homebrew/opt/openssl/include`). +* You can omit the `DGTEST_ROOT` flag to skip tests. + +2\. Make the S2 binary (and associated tests if `GTEST_ROOT` was + specified in CMake: + +
+$ make
+Scanning dependencies of target gtest
+...
+Scanning dependencies of target s2
+...
+[ 35%] Built target s2
+...
+[100%] Built target point_index
+$
+
+ +3\. Run the S2 tests (if `GTEST_ROOT` was specified in CMake): + +
+$ make test
+Running tests...
+Test project /.../s2geometry-master/build
+      Start  1: id_set_lexicon_test
+...
+83/83 Test #83: value_lexicon_test .............................   Passed    0.01 sec
+100% tests passed, 0 tests failed out of 83
+Total Test time (real) =  78.67 sec
+$
+
+ + +4\. Install the built S2 library: + +
+$ make install
+[  1%] Built target gtest
+...
+[100%] Built target point_index
+Install the project...
+-- Install configuration: ""
+...
+$
+
+ +Congratulations! You've successfully installed S2, built the library, and run +all of the tests! You're now ready to move on to the +[C++ Quickstart](/devguide/cpp/quickstart). diff --git a/docs/community/code-of-conduct.md b/docs/community/code-of-conduct.md new file mode 100644 index 00000000..a9ad35e1 --- /dev/null +++ b/docs/community/code-of-conduct.md @@ -0,0 +1,87 @@ +--- +title: Code of Conduct +--- + +This Code of Conduct is adapted from the Contributor Covenant, version 1.4, +available at https://www.contributor-covenant.org/version/1/4/code-of-conduct/ +This Code of Conduct is licensed under the Creative Commons 4.0. + +## Our Pledge + +In the interest of fostering an open and welcoming environment, we as +contributors and maintainers pledge to making participation in our project and +our community a harassment-free experience for everyone, regardless of age, body +size, disability, ethnicity, gender identity and expression, level of +experience, nationality, personal appearance, race, religion, or sexual identity +and orientation. + +## Our Standards + +Examples of behavior that contributes to creating a positive environment +include: + +* Using welcoming and inclusive language +* Being respectful of differing viewpoints and experiences +* Gracefully accepting constructive criticism +* Focusing on what is best for the community +* Showing empathy towards other community members + +Examples of unacceptable behavior by participants include: + +* The use of sexualized language or imagery and unwelcome sexual + attention or advances +* Trolling, insulting/derogatory comments, and personal or political + attacks +* Public or private harassment +* Publishing others’ private information, such as a physical or + electronic address, without explicit permission +* Other conduct which could reasonably be considered inappropriate + in a professional setting + +## Our Responsibilities + +Project maintainers are responsible for clarifying the standards of +acceptable behavior and are expected to take appropriate and fair corrective +action in response to any instances of unacceptable behavior. + +Project maintainers have the right and responsibility to remove, edit, or +reject comments, commits, code, wiki edits, issues, and other contributions +that are not aligned to this Code of Conduct, or to ban temporarily or +permanently any contributor for other behaviors that they deem inappropriate, +threatening, offensive, or harmful. + +## Scope + +This Code of Conduct applies to all individuals associated with the project +community and their actions and behaviors both inside and outside of the project +spaces. + +## Conflict Resolution and Enforcement + +We do not believe that all conflict is bad; healthy debate and disagreement +often yield positive results. However, it is never okay to be disrespectful or +to engage in behavior that violates the project’s code of conduct. + +If you see someone violating the code of conduct, you are encouraged to +address the behavior directly with those involved. Many issues can be resolved +quickly and easily, and this gives people more control over the outcome of their +dispute. If you are unable to resolve the matter for any reason, or if the +behavior is threatening or harassing, report it. We are dedicated to providing +an environment where participants feel welcome and safe. + +Reports should be directed to the +S2 Geometry Library Admin, +the Project Steward for S2. It is the Project Steward’s duty to receive and +address reported violations of the code of conduct. The steward will then work +with a committee consisting of representatives from the Open Source Programs +Office and the Google Cloud Platform Open Source Strategy team. + +We will investigate every complaint, but you may not receive a direct +response. We will use our discretion in determining when and how to follow +up on reported incidents, which may range from not taking action to permanent +expulsion from the project and project-sponsored spaces. We will notify the +accused of the report and provide them an opportunity to discuss it before any +action is taken. In potentially harmful situations, such as ongoing harassment +or threats to anyone's safety, we may take action without notice. The Project +Steward will only affirmatively identify the reporter when it is necessary in +the resolution of the report and with permission from the reporter. diff --git a/docs/community/index.md b/docs/community/index.md new file mode 100644 index 00000000..75f0dd64 --- /dev/null +++ b/docs/community/index.md @@ -0,0 +1,26 @@ +--- +title: Community +--- + +S2 aims to have an active community of developers who are using, enhancing and +building valuable integrations with other software projects. We’d love your help +to improve and extend the project. You can reach us via the S2 Mailing List at +s2geometry-io@googlegroups.com +to start engaging with the project and its members. + +## Code of Conduct + +S2 wants to foster an open and welcoming environment for both our +contributors and maintainers. We pledge to make participation in our project and +our community a harassment-free experience for everyone, regardless of age, body +size, disability, ethnicity, gender identity and expression, level of +experience, nationality, personal appearance, race, religion, or sexual identity +and orientation. To faciliate this open environment, please review our +[Code of Conduct](code-of-conduct). + +## Mailing List + +Any questions or suggestions? Just want to be in the loop of what is going on +with the project? Join the + +S2 Geometry Mailing List. diff --git a/docs/devguide/basic_types.md b/docs/devguide/basic_types.md new file mode 100644 index 00000000..1a6a81d3 --- /dev/null +++ b/docs/devguide/basic_types.md @@ -0,0 +1,1272 @@ +--- +title: Basic Types +--- + + + +## S1Angle + +The `S1Angle` class represents a one-dimensional angle (as opposed to a 2D solid +angle). It has methods for converting angles to or from radians, degrees, and +the E5/E6/E7 representations (i.e. degrees multiplied by 1e5/1e6/1e7 and rounded +to the nearest integer). + +```c++ +class S1Angle { + public: + // These methods construct S1Angle objects from their measure in radians + // or degrees. + static constexpr S1Angle Radians(double radians); + static constexpr S1Angle Degrees(double degrees); + static constexpr S1Angle E5(int32 e5); + static constexpr S1Angle E6(int32 e6); + static constexpr S1Angle E7(int32 e7); + + // The default constructor yields a zero angle. This is useful for STL + // containers and class methods with output arguments. + constexpr S1Angle(); + + // Return an angle larger than any finite angle. + static constexpr S1Angle Infinity(); + + // A explicit shorthand for the default constructor. + static constexpr S1Angle Zero(); + + // Return the angle between two points, which is also equal to the distance + // between these points on the unit sphere. The points do not need to be + // normalized. + S1Angle(const S2Point& x, const S2Point& y); + + // Like the constructor above, but return the angle (i.e., distance) + // between two S2LatLng points. + S1Angle(const S2LatLng& x, const S2LatLng& y); + + constexpr double radians() const; + constexpr double degrees() const; + + int32 e5() const; + int32 e6() const; + int32 e7() const; + + // Return the absolute value of an angle. + S1Angle abs() const; + + // Return the angle normalized to the range (-180, 180] degrees. + S1Angle Normalized() const; + + // Normalize this angle to the range (-180, 180] degrees. + void Normalize(); +}; +``` + +See `s1angle.h` for additional methods, including comparison and arithmetic +operators. For example, if `x` and `y` are angles, then you can write + + if (sin(0.5 * x) > (x + y) / (x - y)) { ... } + +## S2Point + +The `S2Point` class represents a point on the unit sphere as a 3D vector. +Usually points are normalized to be unit length, but some methods do not require +this. The `S2Point` class is simply a synonym for the `Vector_3d` class from +`util/math/vector.h`, which defines overloaded operators that make it convenient +to write arithmetic expressions (e.g. `x*p1 + (1-x)*p2`). + +Some of its more useful methods include: + +```c++ +class S2Point /* Vector3_d */ { + public: + S2Point(double x, double y, double z); + double x(), y(), z(); // Named component accessors + double& operator[](int i); // Return component i (0, 1, 2) + bool operator==(const S2Point& v); // Equality testing + bool operator!=(const S2Point& v); // Inequality testing + S2Point operator+=(const S2Point& v); // Add another vector + S2Point operator-=(const S2Point& v); // Subtract another vector + S2Point operator*=(double k); // Multiply by a scalar + S2Point operator/=(double k); // Divide by a scalar + S2Point operator+(const S2Point& v); // Add two vectors + S2Point operator-(const S2Point& v); // Subtract two vectors + S2Point operator*(double k); // Multiply by a scalar + S2Point operator/(double k); // Divide by a scalar + double DotProd(const S2Point& v); // Dot product + S2Point CrossProd(const S2Point& v); // Cross product + double Norm2(); // Squared L2 norm + double Norm(); // L2 norm + S2Point Normalize(); // Return a normalized **copy** + double Angle(const S2Point&v); // Angle between two vectors (radians) +}; +``` + +Note that the `S2Point` type is technically defined as follows: + + typedef Vector3 Vector3_d; + typedef Vector3_d S2Point; + +See `util/math/vector.h` for details and additional methods. + +## S2LatLng + +The `S2LatLng` class represents a point on the unit sphere as a pair of +latitude-longitude coordinates. + +```c++ +class S2LatLng { + public: + // Constructor. The latitude and longitude are allowed to be outside + // the is_valid() range. However, note that most methods that accept + // S2LatLngs expect them to be normalized (see Normalized() below). + S2LatLng(S1Angle lat, S1Angle lng); + + // The default constructor sets the latitude and longitude to zero. This is + // mainly useful when declaring arrays, STL containers, etc. + S2LatLng(); + + // Convert a direction vector (not necessarily unit length) to an S2LatLng. + explicit S2LatLng(const S2Point& p); + + // Returns an S2LatLng for which is_valid() will return false. + static S2LatLng Invalid(); + + // Convenience functions -- shorter than calling S1Angle::Radians(), etc. + static S2LatLng FromRadians(double lat_radians, double lng_radians); + static S2LatLng FromDegrees(double lat_degrees, double lng_degrees); + static S2LatLng FromE5(int32 lat_e5, int32 lng_e5); + static S2LatLng FromE6(int32 lat_e6, int32 lng_e6); + static S2LatLng FromE7(int32 lat_e7, int32 lng_e7); + + // Accessor methods. + S1Angle lat() const; + S1Angle lng() const; + const R2Point& coords() const; + + // Return true if the latitude is between -90 and 90 degrees inclusive + // and the longitude is between -180 and 180 degrees inclusive. + bool is_valid() const; + + // Clamps the latitude to the range [-90, 90] degrees, and adds or subtracts + // a multiple of 360 degrees to the longitude if necessary to reduce it to + // the range [-180, 180]. + S2LatLng Normalized() const; + + // Convert a normalized S2LatLng to the equivalent unit-length vector. + S2Point ToPoint() const; + + // Return the distance (measured along the surface of the sphere) to the + // given S2LatLng. This is mathematically equivalent to: + // + // S1Angle(ToPoint(), o.ToPoint()) + // + // but this implementation is slightly more efficient. Both S2LatLngs + // must be normalized. + S1Angle GetDistance(const S2LatLng& o) const; +}; +``` + +See `s2latlng.h` for additional methods, including comparison and arithmetic +operators. + +## S2Region + +An `S2Region` represents a two-dimensional region over the unit sphere. It is an +abstract interface with various concrete subtypes, such as discs, rectangles, +polylines, polygons, geometry collections, buffered shapes, etc. + +The main purpose of this interface is to allow complex regions to be +approximated as simpler regions. So rather than having a wide variety of virtual +methods that are implemented by all subtypes, the interface is restricted to +methods that are useful for computing approximations. Here they are: + +```c++ +class S2Region { + public: + virtual ~S2Region(); + + // Return a deep copy of this region. If you want to narrow the result to a + // specific known region type, use down_cast from casts.h. + // Subtypes return pointers to that subtype from their Clone() methods. + virtual S2Region* Clone() const = 0; + + // Return a bounding spherical cap. This is not guaranteed to be exact. + virtual S2Cap GetCapBound() const = 0; + + // Return a bounding latitude-longitude rectangle that contains the region. + // The bounds are not guaranteed to be tight. + virtual S2LatLngRect GetRectBound() const = 0; + + // If this method returns true, the region completely contains the given + // cell. Otherwise, either the region does not contain the cell or the + // containment relationship could not be determined. + virtual bool Contains(const S2Cell& cell) const = 0; + + // If this method returns false, the region does not intersect the given + // cell. Otherwise, either region intersects the cell, or the intersection + // relationship could not be determined. + virtual bool MayIntersect(const S2Cell& cell) const = 0; + + // Return true if and only if the given point is contained by the region. + // The point 'p' is generally required to be unit length, although some + // subtypes may relax this restriction. + virtual bool Contains(const S2Point& p) const = 0; + + ////////////////////////////////////////////////////////////////////////// + // Many S2Region subtypes also define the following non-virtual methods. + + // Appends a serialized representation of the region to "encoder". + void Encode(Encoder* const encoder) const; + + // Decodes an S2Region encoded with Encode(). Note that this method + // requires that an S2Region object of the appropriate concrete type has + // already been constructed. It is not possible to decode regions of + // unknown type. + bool Decode(Decoder* const decoder); +}; +``` + +See `s2region.h` for the full interface. + +## S2LatLngRect + +An `S2LatLngRect` is a type of `S2Region` that represents a rectangle in +latitude-longitude space. It is capable of representing the empty and full +rectangles as well as single points. It has an `AddPoint` method that makes it +easy to construct a bounding rectangle for a set of points, including point sets +that span the 180 degree meridian. Here are its methods: + +```c++ +class S2LatLngRect final : public S2Region { + public: + // Construct a rectangle from minimum and maximum latitudes and longitudes. + // If lo.lng() > hi.lng(), the rectangle spans the 180 degree longitude + // line. Both points must be normalized, with lo.lat() <= hi.lat(). + // The rectangle contains all the points p such that 'lo' <= p <= 'hi', + // where '<=' is defined in the obvious way. + S2LatLngRect(const S2LatLng& lo, const S2LatLng& hi); + + // Construct a rectangle from latitude and longitude intervals. The two + // intervals must either be both empty or both non-empty, and the latitude + // interval must not extend outside [-90, +90] degrees. + S2LatLngRect(const R1Interval& lat, const S1Interval& lng); + + // The default constructor creates an empty S2LatLngRect. + S2LatLngRect(); + + // Construct a rectangle of the given size centered around the given point. + static S2LatLngRect FromCenterSize(const S2LatLng& center, + const S2LatLng& size); + + // Construct a rectangle containing a single (normalized) point. + static S2LatLngRect FromPoint(const S2LatLng& p); + + // Construct the minimal bounding rectangle containing the two given + // normalized points. The points do not need to be ordered. + static S2LatLngRect FromPointPair(const S2LatLng& p1, const S2LatLng& p2); + + // Accessor methods. + S1Angle lat_lo() const; + S1Angle lat_hi() const; + S1Angle lng_lo() const; + S1Angle lng_hi() const; + const R1Interval& lat() const; + const S1Interval& lng() const; + R1Interval *mutable_lat(); + S1Interval *mutable_lng(); + S2LatLng lo() const; + S2LatLng hi() const; + + // The canonical empty and full rectangles. + static S2LatLngRect Empty(); + static S2LatLngRect Full(); + + // Return true if the rectangle is valid, which essentially just means + // that the latitude bounds do not exceed Pi/2 in absolute value and + // the longitude bounds do not exceed Pi in absolute value. Also, if + // either the latitude or longitude bound is empty then both must be. + bool is_valid() const; + + // Return true if the rectangle is empty, i.e. it contains no points at all. + bool is_empty() const; + + // Return true if the rectangle is full, i.e. it contains all points. + bool is_full() const; + + // Return true if the rectangle is a point, i.e. lo() == hi() + bool is_point() const; + + // Return true if the rectangle crosses the 180 degree longitude line. + bool is_inverted() const; + + // Return the k-th vertex of the rectangle (k = 0,1,2,3) in CCW order (lower + // left, lower right, upper right, upper left). + S2LatLng GetVertex(int k) const; + + // Return the center of the rectangle in latitude-longitude space + // (in general this is not the center of the region on the sphere). + S2LatLng GetCenter() const; + + // Return the width and height of this rectangle in latitude-longitude + // space. Empty rectangles have a negative width and height. + S2LatLng GetSize() const; + + // Returns the surface area of this rectangle on the unit sphere. + double Area() const; + + // Return the true centroid of the rectangle multiplied by its surface area. + // The true centroid is the "center of mass" of the rectangle viewed as + // subset of the unit sphere, i.e. it is the point in space about which this + // curved shape would rotate.) + // + // The reason for multiplying the result by the rectangle area is to make it + // easier to compute the centroid of more complicated shapes. + S2Point GetCentroid() const; + + // More efficient version of Contains() that accepts a S2LatLng rather than + // an S2Point. + bool Contains(const S2LatLng& ll) const; + + // Return true if and only if the given point is contained in the interior + // of the region (i.e. the region excluding its boundary). + bool InteriorContains(const S2Point& p) const; + + // More efficient version of InteriorContains() that accepts a S2LatLng + // rather than an S2Point. + bool InteriorContains(const S2LatLng& ll) const; + + // Return true if and only if the rectangle contains the given other + // rectangle. + bool Contains(const S2LatLngRect& other) const; + + // Return true if and only if the interior of this rectangle contains all + // points of the given other rectangle (including its boundary). + bool InteriorContains(const S2LatLngRect& other) const; + + // Return true if this rectangle and the given other rectangle have any + // points in common. + bool Intersects(const S2LatLngRect& other) const; + + // Returns true if this rectangle intersects the given cell. (This is an + // exact test and may be fairly expensive, see also MayIntersect below.) + bool Intersects(const S2Cell& cell) const; + + // Return true if and only if the interior of this rectangle intersects + // any point (including the boundary) of the given other rectangle. + bool InteriorIntersects(const S2LatLngRect& other) const; + + // Increase the size of the bounding rectangle to include the given point. + // The rectangle is expanded by the minimum amount possible. + void AddPoint(const S2Point& p); + void AddPoint(const S2LatLng& ll); + + // Return a rectangle that has been expanded by margin.lat() on each side in + // the latitude direction, and by margin.lng() on each side in the longitude + // direction. If either margin is negative, then shrink the rectangle on + // the corresponding sides instead. The resulting rectangle may be empty. + // + // If you are trying to grow a rectangle by a certain *distance* on the + // sphere (e.g. 5km), use the ExpandedByDistance() method instead. + S2LatLngRect Expanded(const S2LatLng& margin) const; + + // If the rectangle does not include either pole, return it unmodified. + // Otherwise expand the longitude range to Full() so that the rectangle + // contains all possible representations of the contained pole(s). + S2LatLngRect PolarClosure() const; + + // Return the smallest rectangle containing the union of this rectangle and + // the given rectangle. + S2LatLngRect Union(const S2LatLngRect& other) const; + + // Return the smallest rectangle containing the intersection of this + // rectangle and the given rectangle. Note that the region of intersection + // may consist of two disjoint rectangles, in which case a single rectangle + // spanning both of them is returned. + S2LatLngRect Intersection(const S2LatLngRect& other) const; + + // Expand this rectangle so that it contains all points within the given + // distance of the boundary, and return the smallest such rectangle. If the + // distance is negative, then instead shrink this rectangle so that it + // excludes all points within the given absolute distance of the boundary, + // and return the largest such rectangle. + // + // Unlike Expanded(), this method treats the rectangle as a set of points on + // the sphere, and measures distances on the sphere. For example, you can + // use this method to find a rectangle that contains all points within 5km + // of a given rectangle. + S2LatLngRect ExpandedByDistance(S1Angle distance) const; + + // Returns the minimum distance (measured along the surface of the sphere) to + // the given S2LatLngRect. Both S2LatLngRects must be non-empty. + S1Angle GetDistance(const S2LatLngRect& other) const; + + // Returns the minimum distance (measured along the surface of the sphere) + // from a given point to the rectangle (both its boundary and its interior). + S1Angle GetDistance(const S2LatLng& p) const; + + //////////////////////////////////////////////////////////////////////// + // S2Region interface (see s2region.h for details): + + S2LatLngRect* Clone() const override; + S2Cap GetCapBound() const override; + S2LatLngRect GetRectBound() const override; + // The point 'p' does not need to be normalized. + bool Contains(const S2Cell& cell) const override; + + // This test is cheap but is NOT exact. Use Intersects() if you want a more + // accurate and more expensive test. + bool MayIntersect(const S2Cell& cell) const override; +}; +``` + +See `s2latlngrect.h` for additional methods. + +## S1ChordAngle + +The S2 library actually defines two angle representations, `S1Angle` and +`S1ChordAngle`. + +`S1ChordAngle` is a specialized angle type that represents the distance +between two points on the sphere. (Note that the distance between two +points can be expressed as the angle between those points measured from the +sphere's center.) Its representation makes it very efficient for computing +and comparing distances, but unlike `S1Angle` it is only capable of +representing angles between 0 and 180 degrees. (Internally, `S1ChordAngle` +computes the squared distance between two points through the interior of +the sphere, i.e. the squared length of the chord between those points). + +`S1ChordAngle` is the preferred representation for distances internally, +because it is much faster to compute than `S1Angle` and also supports the +exact predicates needed for robust geometric algorithms. `S1ChordAngle` +supports many of the same methods as `S1Angle`, and it is also easy to +convert back and forth as required. + +## S2Cap + +An `S2Cap` represents a spherical cap, i.e. a portion of a sphere cut off by a +plane. This is the equivalent of a *closed disc* in planar geometry (i.e., a +circle together with its interior), so if you are looking for a way to represent +a circle or disc then you are probably looking for an `S2Cap`. + +Here are the methods available: + +```c++ +class S2Cap final : public S2Region { + public: + // The default constructor returns an empty S2Cap. + S2Cap(); + + // Construct a cap with the given center and radius. A negative radius + // yields an empty cap; a radius of 180 degrees or more yields a full cap + // (containing the entire sphere). + S2Cap(const S2Point& center, S1Angle radius); + + // Convenience function that creates a cap containing a single point. This + // method is more efficient that the S2Cap(center, radius) constructor. + static S2Cap FromPoint(const S2Point& center); + + // Return a cap with the given center and height. (The height of a cap is + // the distance from the cap center to the cutoff plane; see comments + // above.) A negative height yields an empty cap; a height of 2 or more + // yields a full cap. + static S2Cap FromCenterHeight(const S2Point& center, double height); + + // Return a cap with the given center and surface area. Note that the area + // can also be interpreted as the solid angle subtended by the cap (because + // the sphere has unit radius). A negative area yields an empty cap; an + // area of 4*Pi or more yields a full cap. + static S2Cap FromCenterArea(const S2Point& center, double area); + + // Return an empty cap, i.e. a cap that contains no points. + static S2Cap Empty(); + + // Return a full cap, i.e. a cap that contains all points. + static S2Cap Full(); + + // Accessor methods. + const S2Point& center() const; + S1ChordAngle radius() const; + double height() const; + + // Return the cap radius. (This method is relatively expensive since only + // the cap height is stored internally, and may yield a slightly different + // result than the value passed to the (center, radius) constructor.) + S1Angle GetRadius() const; + + // Return the area of the cap. + double GetArea() const; + + // Return the true centroid of the cap multiplied by its surface area. Note + // that if you just want the "surface centroid", i.e. the normalized result, + // then it is much simpler just to call center(). + // + // The reason for multiplying the result by the cap area is to make it + // easier to compute the centroid of more complicated shapes. + S2Point GetCentroid() const; + + // We allow negative heights (to represent empty caps) but heights are + // normalized so that they do not exceed 2. + bool is_valid() const; + + // Return true if the cap is empty, i.e. it contains no points. + bool is_empty() const; + + // Return true if the cap is full, i.e. it contains all points. + bool is_full() const; + + // Return the complement of the interior of the cap. A cap and its + // complement have the same boundary but do not share any interior points. + S2Cap Complement() const; + + // Return true if and only if this cap contains the given other cap + // (in a set containment sense, e.g. every cap contains the empty cap). + bool Contains(const S2Cap& other) const; + + // Return true if and only if this cap intersects the given other cap, + // i.e. whether they have any points in common. + bool Intersects(const S2Cap& other) const; + + // Return true if and only if the interior of this cap intersects the + // given other cap. (This relationship is not symmetric, since only + // the interior of this cap is used.) + bool InteriorIntersects(const S2Cap& other) const; + + // Return true if and only if the given point is contained in the interior + // of the cap (i.e. the cap excluding its boundary). + bool InteriorContains(const S2Point& p) const; + + // Increase the cap height if necessary to include the given point. If the + // cap is empty then the center is set to the given point, but otherwise the + // center is not changed. + void AddPoint(const S2Point& p); + + // Increase the cap height if necessary to include "other". If the current + // cap is empty it is set to the given other cap. + void AddCap(const S2Cap& other); + + // Return a cap that contains all points within a given distance of this + // cap. Note that any expansion of the empty cap is still empty. + S2Cap Expanded(S1Angle distance) const; + + // Return the smallest cap which encloses this cap and "other". + S2Cap Union(const S2Cap& other) const; + + //////////////////////////////////////////////////////////////////////// + // S2Region interface (see s2region.h for details): + + S2Cap* Clone() const override; + S2Cap GetCapBound() const override; + S2LatLngRect GetRectBound() const override; + // The point "p" should be a unit-length vector. + bool Contains(const S2Cell& cell) const override; + bool MayIntersect(const S2Cell& cell) const override; +}; +``` + +See `s2cap.h` for additional methods. + +## S2Polyline + +An `S2Polyline` represents a sequence of zero or more vertices connected by +straight edges (geodesics). Edges of length 0 and 180 degrees are not allowed, +i.e. adjacent vertices should not be identical or antipodal. + +```c++ +class S2Polyline final : public S2Region { + public: + // Creates an empty S2Polyline that should be initialized by calling Init() + // or Decode(). + S2Polyline(); + + // Convenience constructors that call Init() with the given vertices. + explicit S2Polyline(const std::vector& vertices); + explicit S2Polyline(const std::vector& vertices); + + ~S2Polyline(); + + // Initialize a polyline that connects the given vertices. Empty polylines are + // allowed. Adjacent vertices should not be identical or antipodal. All + // vertices should be unit length. + void Init(const std::vector& vertices); + + // Convenience initialization function that accepts latitude-longitude + // coordinates rather than S2Points. + void Init(const std::vector& vertices); + + // Return true if the given vertices form a valid polyline. + bool IsValid() const; + + // Returns true if this is *not* a valid polyline and sets "error" + // appropriately. Otherwise returns false and leaves "error" unchanged. + bool FindValidationError(S2Error* error) const; + + // Accessor methods. + int num_vertices() const; + const S2Point& vertex(int k) const; + + // Return the length of the polyline. + S1Angle GetLength() const; + + // Return the true centroid of the polyline multiplied by the length of the + // polyline. Prescaling by the polyline length makes it easy to compute the + // centroid of several polylines (by simply adding up their centroids). + S2Point GetCentroid() const; + + // Return the point whose distance from vertex 0 along the polyline is the + // given fraction of the polyline's total length. Fractions less than zero + // or greater than one are clamped. The return value is unit length. + S2Point Interpolate(double fraction) const; + + // Like Interpolate(), but also return the index of the next polyline + // vertex after the interpolated point P. This allows the caller to easily + // construct a given suffix of the polyline by concatenating P with the + // polyline vertices starting at "next_vertex". + S2Point GetSuffix(double fraction, int* next_vertex) const; + + // The inverse operation of GetSuffix/Interpolate. Given a point on the + // polyline, returns the ratio of the distance to the point from the + // beginning of the polyline over the length of the polyline. The return + // value is always betwen 0 and 1 inclusive. See GetSuffix() for the + // meaning of "next_vertex". + double UnInterpolate(const S2Point& point, int next_vertex) const; + + // Given a point, returns a point on the polyline that is closest to the given + // point. See GetSuffix() for the meaning of "next_vertex". + S2Point Project(const S2Point& point, int* next_vertex) const; + + // Returns true if the point given is on the right hand side of the polyline, + // using a naive definition of "right-hand-sideness" where the point is on + // the RHS of the polyline iff the point is on the RHS of the line segment in + // the polyline which it is closest to. + bool IsOnRight(const S2Point& point) const; + + // Return true if this polyline intersects the given polyline. If the + // polylines share a vertex they are considered to be intersecting. + bool Intersects(const S2Polyline* line) const; + + // Reverse the order of the polyline vertices. + void Reverse(); + + // Return a subsequence of vertex indices such that the polyline connecting + // these vertices is never further than "tolerance" from the original + // polyline. + void SubsampleVertices(S1Angle tolerance, std::vector* indices) const; + + // Return true if "covered" is within "max_error" of a contiguous subpath of + // this polyline over its entire length. Specifically, this method returns + // true if this polyline has parameterization a:[0,1] -> S^2, "covered" has + // parameterization b:[0,1] -> S^2, and there is a non-decreasing function + // f:[0,1] -> [0,1] such that distance(a(f(t)), b(t)) <= max_error for all t. + // + // You can think of this as testing whether it is possible to drive a car + // along "covered" and a car along some subpath of this polyline such that no + // car ever goes backward, and the cars are always within "max_error" of each + // other. + bool NearlyCoversPolyline(const S2Polyline& covered, + S1Angle max_error) const; + + //////////////////////////////////////////////////////////////////////// + // S2Region interface (see s2region.h for details): + + S2Polyline* Clone() const override; + S2Cap GetCapBound() const override; + S2LatLngRect GetRectBound() const override; + bool Contains(const S2Cell& cell) const override; + bool MayIntersect(const S2Cell& cell) const override; +}; +``` + +## S2Loop + +An `S2Loop` represents a simple spherical polygon. It consists of a single chain +of vertices where the first vertex is implicitly connected to the last. All +loops are defined to have a CCW orientation, i.e. the interior of the loop is on +the left side of the edges. This implies that a clockwise loop enclosing a small +area is interpreted to be a CCW loop enclosing a very large area. + +Loops are not allowed to have any duplicate vertices (whether adjacent or not), +and non-adjacent edges are not allowed to intersect. Loops must have at least 3 +vertices (except for the "empty" and "full" loops discussed below). These +restrictions make it possible to implement exact polygon-polygon containment and +intersection tests very efficiently. See `S2Builder` if your data does not meet +these requirements. + +There are two special loops: the "empty" loop contains no points, while the +"full" loop contains all points. These loops do not have any edges, but to +preserve the invariant that every loop can be represented as a vertex chain, +they are defined as having exactly one vertex each (see kEmpty and kFull). + +Point containment of loops is defined such that if the sphere is subdivided into +faces (loops), every point is contained by exactly one face. This implies that +loops do not necessarily contain their vertices. + +```c++ +class S2Loop final : public S2Region { + public: + // Default constructor. The loop must be initialized by calling Init() or + // Decode() before it is used. + S2Loop(); + + // Convenience constructor that calls Init() with the given vertices. + explicit S2Loop(const std::vector& vertices); + + // Construct a loop corresponding to the given cell. + // + // Note that the loop and cell *do not* contain exactly the same set of + // points, because S2Loop and S2Cell have slightly different definitions of + // point containment. For example, an S2Cell vertex is contained by all + // four neighboring S2Cells, but it is contained by exactly one of four + // S2Loops constructed from those cells. + explicit S2Loop(const S2Cell& cell); + + ~S2Loop(); + + // Initialize a loop with given vertices. The last vertex is implicitly + // connected to the first. All points should be unit length. Loops must + // have at least 3 vertices (except for the "empty" and "full" loops, see + // kEmpty and kFull). This method may be called multiple times. + void Init(const std::vector& vertices); + + // A special vertex chain of length 1 that creates an empty loop (i.e., a + // loop with no edges that contains no points). Example usage: + // S2Loop empty(S2Loop::kEmpty()); + static std::vector kEmpty(); + + // A special vertex chain of length 1 that creates a full loop (i.e., a loop + // with no edges that contains all points). See kEmpty() for details. + static std::vector kFull(); + + // Returns true if this is a valid loop. + bool IsValid() const; + + // Returns true if this is *not* a valid loop and sets "error" + // appropriately. Otherwise returns false and leaves "error" unchanged. + bool FindValidationError(S2Error* error) const; + + int num_vertices() const; + + // For convenience, we make two entire copies of the vertex list available: + // vertex(n..2*n-1) is mapped to vertex(0..n-1), where n == num_vertices(). + const S2Point& vertex(int i) const; + + // Like vertex(), but this method returns vertices in reverse order if the + // loop represents a polygon hole. This ensures that the interior of the + // polygon is always to the left of the vertex chain. + const S2Point& oriented_vertex(int i) const; + + // Return true if this is the special "empty" loop that contains no points. + bool is_empty() const; + + // Return true if this is the special "full" loop that contains all points. + bool is_full() const; + + // Return true if this loop is either "empty" or "full". + bool is_empty_or_full() const; + + // The depth of a loop is defined as its nesting level within its containing + // polygon. "Outer shell" loops have depth 0, holes within those loops have + // depth 1, shells within those holes have depth 2, etc. This field is only + // used by the S2Polygon implementation. + int depth() const; + void set_depth(int depth); + + // Return true if this loop represents a hole in its containing polygon. + bool is_hole() const; + + // The sign of a loop is -1 if the loop represents a hole in its containing + // polygon, and +1 otherwise. + int sign() const; + + // Return true if the loop area is at most 2*Pi. + bool IsNormalized() const; + + // Invert the loop if necessary so that the area enclosed by the loop is at + // most 2*Pi. + void Normalize(); + + // Reverse the order of the loop vertices, effectively complementing + // the region represented by the loop. + void Invert(); + + // Return the area of the loop interior, i.e. the region on the left side of + // the loop. The return value is between 0 and 4*Pi. + double GetArea() const; + + // Return the true centroid of the loop multiplied by the area of the loop. + // We prescale by the loop area for two reasons: (1) it is cheaper to + // compute this way, and (2) it makes it easier to compute the centroid of + // more complicated shapes (by splitting them into disjoint regions and + // adding their centroids). + S2Point GetCentroid() const; + + // Return the sum of the turning angles at each vertex. The return value is + // positive if the loop is counter-clockwise, negative if the loop is + // clockwise, and zero if the loop is a great circle. + // + // This quantity is also called the "geodesic curvature" of the loop. + double GetCurvature() const; + + // Return the maximum error in GetCurvature(). The return value is not + // constant; it depends on the loop. + double GetCurvatureMaxError() const; + + // Return the distance from the given point to the loop interior. If the + // loop is empty, return S1Angle::Infinity(). + S1Angle GetDistance(const S2Point& x) const; + + // Return the distance from the given point to the loop boundary. If the + // loop is empty or full, return S1Angle::Infinity() (since the loop has no + // boundary). + S1Angle GetDistanceToBoundary(const S2Point& x) const; + + // If the given point is contained by the loop, return it. Otherwise return + // the closest point on the loop boundary. If the loop is empty, return the + // input argument. + S2Point Project(const S2Point& x) const; + + // Return the closest point on the loop boundary to the given point. If the + // loop is empty or full, return the input argument (since the loop has no + // boundary). + S2Point ProjectToBoundary(const S2Point& x) const; + + // Return true if the region contained by this loop is a superset of the + // region contained by the given other loop. + bool Contains(const S2Loop* b) const; + + // Return true if the region contained by this loop intersects the region + // contained by the given other loop. + bool Intersects(const S2Loop* b) const; + + // Return true if two loops have the same vertices in the same linear order + // (i.e., cyclic rotations are not allowed). + bool Equals(const S2Loop* b) const; + + // Return true if two loops have the same boundary. This is true if and + // only if the loops have the same vertices in the same cyclic order (i.e., + // the vertices may be cyclically rotated). The empty and full loops are + // considered to have different boundaries. + bool BoundaryEquals(const S2Loop* b) const; + + // Return true if two loops have the same boundary except for vertex + // perturbations. More precisely, the vertices in the two loops must be in + // the same cyclic order, and corresponding vertex pairs must be separated + // by no more than "max_error". + bool BoundaryApproxEquals(const S2Loop& b, + S1Angle max_error = S1Angle::Radians(1e-15)) const; + + // Return true if the two loop boundaries are within "max_error" of each + // other along their entire lengths. The two loops may have different + // numbers of vertices. More precisely, this method returns true if the two + // loops have parameterizations a:[0,1] -> S^2, b:[0,1] -> S^2 such that + // distance(a(t), b(t)) <= max_error for all t. You can think of this as + // testing whether it is possible to drive two cars all the way around the + // two loops such that no car ever goes backward and the cars are always + // within "max_error" of each other. + bool BoundaryNear(const S2Loop& b, + S1Angle max_error = S1Angle::Radians(1e-15)) const; + + // This method computes the oriented surface integral of some quantity f(x) + // over the loop interior, given a function f_tri(A,B,C) that returns the + // corresponding integral over the spherical triangle ABC. + template + T GetSurfaceIntegral(T f_tri(const S2Point&, const S2Point&, const S2Point&)) + const; + + //////////////////////////////////////////////////////////////////////// + // S2Region interface (see s2region.h for details): + + S2Loop* Clone() const override; + + // GetRectBound() returns essentially tight results, while GetCapBound() + // might have a lot of extra padding. Both bounds are conservative in that + // if the loop contains a point P, then the bound contains P also. + S2Cap GetCapBound() const override; + S2LatLngRect GetRectBound() const override; + bool Contains(const S2Cell& cell) const override; + bool MayIntersect(const S2Cell& cell) const override; + // Return true if the loop contains the given point. Point containment is + // defined such that if the sphere is subdivided into faces (loops), every + // point is contained by exactly one face. This implies that loops do not + // necessarily contain their vertices. + bool Contains(const S2Point& p) const override; + + // Generally clients should not use S2Loop::Encode(). Instead they should + // encode an S2Polygon, which unlike this method supports (lossless) + // compression. + void Encode(Encoder* const encoder) const override; + + // Decode a loop encoded with Encode() or EncodeCompressed(). These methods + // may be called with loops that have already been initialized. + bool Decode(Decoder* const decoder) override; + bool DecodeWithinScope(Decoder* const decoder) override; +}; +``` + +## S2Polygon + +An `S2Polygon` is an `S2Region` object that represents a polygon. A polygon is +defined by zero or more loops; recall that the interior of a loop is defined to +be its left-hand side (see `S2Loop`). There are two different conventions for +creating an `S2Polygon`: + +* `InitNested()` expects the input loops to be nested hierarchically. The + polygon interior then consists of the set of points contained by an odd + number of loops. So for example, a circular region with a hole in it would + be defined as two CCW loops, with one loop containing the other. The loops + can be provided in any order. + + When the orientation of the input loops is unknown, the nesting requirement + is typically met by calling `S2Loop::Normalize()` on each loop (which + inverts the loop if necessary so that it encloses at most half the sphere). + But in fact any set of loops can be used as long as (1) there is no pair of + loops that cross, and (2) there is no pair of loops whose union is the + entire sphere. + +* `InitOriented()` expects the input loops to be oriented such that the + polygon interior is on the left-hand side of every loop. So for example, a + circular region with a hole in it would be defined using a CCW outer loop + and a CW inner loop. The loop orientations must all be consistent; for + example, it is not valid to have one CCW loop nested inside another CCW + loop, because the region between the two loops is on the left-hand side of + one loop and the right-hand side of the other. + +Most clients will not call these methods directly; instead they should use +`S2Builder`, which has better support for dealing with imperfect data. + +When the polygon is initialized, the given loops are automatically converted +into a canonical form consisting of "shells" and "holes". Shells and holes are +both oriented CCW, and are nested hierarchically. The loops are reordered to +correspond to a preorder traversal of the nesting hierarchy; `InitOriented()` +may also invert some loops. + +Polygons may represent any region of the sphere with a polygonal boundary, +including the entire sphere (known as the "full" polygon). The full polygon +consists of a single full loop (see `S2Loop`), whereas the empty polygon has no +loops at all. + +Polygons have the following restrictions: + +* Loops may not cross, i.e. the boundary of a loop may not intersect both the + interior and exterior of any other loop. + +* Loops may not share edges, i.e. if a loop contains an edge AB, then no other + loop may contain AB or BA. + +* Loops may share vertices, however no vertex may appear twice in a single + loop (see `S2Loop`). + +* No loop may be empty. The full loop may appear only in the full polygon. + + +```c++ +class S2Polygon final : public S2Region { + public: + // The default constructor creates an empty polygon. It can be made + // non-empty by calling Init(), Decode(), etc. + S2Polygon(); + + // Convenience constructor that calls InitNested() with the given loops. + // Takes ownership of the loops and clears the vector. + explicit S2Polygon(std::vector> loops); + + // Convenience constructor that creates a polygon with a single loop + // corresponding to the given cell. + explicit S2Polygon(const S2Cell& cell); + + // Convenience constructor that calls Init(unique_ptr). + explicit S2Polygon(std::unique_ptr loop); + + // Destroys the polygon and frees its loops. + ~S2Polygon(); + + // Create a polygon from a set of hierarchically nested loops. The polygon + // interior consists of the points contained by an odd number of loops. + // (Recall that a loop contains the set of points on its left-hand side.) + // + // This method takes ownership of the given loops and clears the given + // vector. It then figures out the loop nesting hierarchy and assigns every + // loop a depth. Shells have even depths, and holes have odd depths. Note + // that the loops are reordered so the hierarchy can be traversed more + // easily (see GetParent(), GetLastDescendant(), and S2Loop::depth()). + // + // This method may be called more than once, in which case any existing + // loops are deleted before being replaced by the input loops. + void InitNested(std::vector> loops); + + // Like InitNested(), but expects loops to be oriented such that the polygon + // interior is on the left-hand side of all loops. This implies that shells + // and holes should have opposite orientations in the input to this method. + void InitOriented(std::vector> loops); + + // Initialize a polygon from a single loop. Takes ownership of the loop. + void Init(std::unique_ptr loop); + + // Releases ownership of and returns the loops of this polygon, and resets + // the polygon to be empty. + std::vector> Release(); + + // Makes a deep copy of the given source polygon. The destination polygon + // will be cleared if necessary. + void Copy(const S2Polygon* src); + + // Returns true if this is a valid polygon (including checking whether all + // the loops are themselves valid). + bool IsValid() const; + + // Returns true if this is *not* a valid polygon and sets "error" + // appropriately. Otherwise returns false and leaves "error" unchanged. + bool FindValidationError(S2Error* error) const; + + // Return true if this is the empty polygon (consisting of no loops). + bool is_empty() const; + + // Return true if this is the full polygon (consisting of a single loop that + // encompasses the entire sphere). + bool is_full() const; + + // Return the number of loops in this polygon. + int num_loops() const; + + // Total number of vertices in all loops. + int num_vertices() const; + + // Return the loop at the given index. Note that during initialization, the + // given loops are reordered according to a preorder traversal of the loop + // nesting hierarchy. This implies that every loop is immediately followed + // by its descendants. This hierarchy can be traversed using the methods + // GetParent(), GetLastDescendant(), and S2Loop::depth(). + const S2Loop* loop(int k) const; + S2Loop* loop(int k); + + // Return the index of the parent of loop k, or -1 if it has no parent. + int GetParent(int k) const; + + // Return the index of the last loop that is contained within loop k. + int GetLastDescendant(int k) const; + + // Return the area of the polygon interior, i.e. the region on the left side + // of an odd number of loops. The return value is between 0 and 4*Pi. + double GetArea() const; + + // Return the true centroid of the polygon multiplied by the area of the + // polygon. We prescale by the polygon area for two reasons: (1) it is + // cheaper to compute this way, and (2) it makes it easier to compute the + // centroid of more complicated shapes (by splitting them into disjoint + // regions and adding their centroids). + S2Point GetCentroid() const; + + // If all of the polygon's vertices happen to be the centers of S2Cells at + // some level, then return that level, otherwise return -1. + int GetSnapLevel() const; + + // Return the distance from the given point to the polygon interior. If the + // polygon is empty, return S1Angle::Infinity(). + S1Angle GetDistance(const S2Point& x) const; + + // Return the distance from the given point to the polygon boundary. If the + // polygon is empty or full, return S1Angle::Infinity() (since the polygon + // has no boundary). + S1Angle GetDistanceToBoundary(const S2Point& x) const; + + // Return the overlap fractions between two polygons, i.e. the ratios of the + // area of intersection to the area of each polygon. + static std::pair GetOverlapFractions(const S2Polygon* a, + const S2Polygon* b); + + // If the given point is contained by the polygon, return it. Otherwise + // return the closest point on the polygon boundary. If the polygon is + // empty, return the input argument. + S2Point Project(const S2Point& x) const; + + // Return the closest point on the polygon boundary to the given point. If + // the polygon is empty or full, return the input argument (since the + // polygon has no boundary). + S2Point ProjectToBoundary(const S2Point& x) const; + + // Return true if this polygon contains the given other polygon, i.e. + // if polygon A contains all points contained by polygon B. + bool Contains(const S2Polygon* b) const; + + // Returns true if this polgyon (A) approximately contains the given other + // polygon (B). This is true if it is possible to move the vertices of B + // no further than "tolerance" such that A contains the modified B. + // + // For example, the empty polygon will contain any polygon whose maximum + // width is no more than "tolerance". + bool ApproxContains(const S2Polygon* b, S1Angle tolerance) const; + + // Return true if this polygon intersects the given other polygon, i.e. + // if there is a point that is contained by both polygons. + bool Intersects(const S2Polygon* b) const; + + // Returns true if this polgyon (A) and the given polygon (B) are + // approximately disjoint. This is true if it is possible to ensure that A + // and B do not intersect by moving their vertices no further than + // "tolerance". + // + // For example, any polygon is approximately disjoint from a polygon whose + // maximum width is no more than "tolerance". + bool ApproxDisjoint(const S2Polygon& b, S1Angle tolerance) const; + + // Initialize this polygon to the intersection, union, difference (A - B), + // or symmetric difference (XOR) of the given two polygons. + // + // "snap_function" allows you to specify a minimum spacing between output + // vertices, and/or that the vertices should be snapped to a discrete set of + // points (e.g. S2CellId centers or E7 lat/lng coordinates). Any snap + // function can be used, including the IdentitySnapFunction with a + // snap_radius of zero (which preserves the input vertices exactly). + void InitToIntersection(const S2Polygon* a, const S2Polygon* b); + void InitToIntersection(const S2Polygon& a, const S2Polygon& b, + const S2Builder::SnapFunction& snap_function); + + void InitToUnion(const S2Polygon* a, const S2Polygon* b); + void InitToUnion(const S2Polygon& a, const S2Polygon& b, + const S2Builder::SnapFunction& snap_function); + + void InitToDifference(const S2Polygon* a, const S2Polygon* b); + void InitToDifference(const S2Polygon& a, const S2Polygon& b, + const S2Builder::SnapFunction& snap_function); + + void InitToSymmetricDifference(const S2Polygon* a, const S2Polygon* b); + void InitToSymmetricDifference( + const S2Polygon& a, const S2Polygon& b, + const S2Builder::SnapFunction& snap_function); + + // Convenience functions that use the IdentitySnapFunction with the given + // snap radius. + void InitToApproxIntersection(const S2Polygon* a, const S2Polygon* b, + S1Angle vertex_merge_radius); + void InitToApproxUnion(const S2Polygon* a, const S2Polygon* b, + S1Angle vertex_merge_radius); + void InitToApproxDifference(const S2Polygon* a, const S2Polygon* b, + S1Angle vertex_merge_radius); + void InitToApproxSymmetricDifference(const S2Polygon* a, const S2Polygon* b, + S1Angle vertex_merge_radius); + + // Initializes this polygon according to the given "snap_function" and + // reduces the number of vertices if possible, while ensuring that no vertex + // moves further than snap_function.snap_radius(). + // + // Simplification works by replacing nearly straight chains of short edges + // with longer edges, in a way that preserves the topology of the input + // polygon up to the creation of degeneracies. This means that loops or + // portions of loops may become degenerate, in which case they are removed. + void InitToSimplified(const S2Polygon& a, + const S2Builder::SnapFunction& snap_function); + + // Use S2Builder to build this polygon by assembling the edges of a + // given polygon after snapping its vertices to the center of leaf cells at + // the given "snap_level". The default snap level corresponds to a + // tolerance of approximately 1.5cm on the surface of the Earth. + // Such a polygon can be efficiently compressed when serialized. + void InitToSnapped(const S2Polygon* polygon, + int snap_level = S2CellId::kMaxLevel); + + // Initialize this polygon to the complement of the given polygon. + void InitToComplement(const S2Polygon* a); + + // Invert the polygon (replace it by its complement). + void Invert(); + + // Intersect this polygon with the polyline "in" and append the resulting + // zero or more polylines to "out" (which must be empty). The polylines + // are appended in the order they would be encountered by traversing "in" + // from beginning to end. Note that the output may include polylines with + // only one vertex, but there will not be any zero-vertex polylines. + std::vector> + IntersectWithPolyline(const S2Polyline& in) const; + + // Similar to IntersectWithPolyline(), except that vertices will be + // dropped as necessary to ensure that all adjacent vertices in the + // sequence obtained by concatenating the output polylines will be + // farther than "vertex_merge_radius" apart. + std::vector> + ApproxIntersectWithPolyline(const S2Polyline& in, + S1Angle vertex_merge_radius) const; + + // Like IntersectWithPolyline, but subtracts this polygon from + // the given polyline. + std::vector> + SubtractFromPolyline(const S2Polyline& in) const; + + // Like ApproxIntersectWithPolyline, but subtracts this polygon + // from the given polyline. + std::vector> + ApproxSubtractFromPolyline(const S2Polyline& in, + S1Angle vertex_merge_radius) const; + + // Return a polygon which is the union of the given polygons. + // Clears the vector and deletes the polygons. + static std::unique_ptr + DestructiveUnion(std::vector> polygons); + static std::unique_ptr + DestructiveApproxUnion(std::vector> polygons, + S1Angle vertex_merge_radius); + + // Initialize this polygon to the outline of the given cell union. + // In principle this polygon should exactly contain the cell union and + // this polygon's inverse should not intersect the cell union, but rounding + // issues may cause this not to be the case. + void InitToCellUnionBorder(const S2CellUnion& cells); + + // Return true if every loop of this polygon shares at most one vertex with + // its parent loop. Every polygon has a unique normalized form. Normalized + // polygons are useful for testing since it is easy to compare whether two + // polygons represent the same region. + bool IsNormalized() const; + + // Return true if two polygons have exactly the same loops. The loops must + // appear in the same order, and corresponding loops must have the same + // linear vertex ordering (i.e., cyclic rotations are not allowed). + bool Equals(const S2Polygon* b) const; + + // Return true if two polygons have the same boundary. More precisely, this + // method requires that both polygons have loops with the same cyclic vertex + // order and the same nesting hierarchy. (This implies that vertices may be + // cyclically rotated between corresponding loops, and the loop ordering may + // be different between the two polygons as long as the nesting hierarchy is + // the same.) + bool BoundaryEquals(const S2Polygon* b) const; + + // Return true if two polygons have the same boundary except for vertex + // perturbations. Both polygons must have loops with the same cyclic vertex + // order and the same nesting hierarchy, but the vertex locations are + // allowed to differ by up to "max_error". + bool BoundaryApproxEquals(const S2Polygon& b, + S1Angle max_error = S1Angle::Radians(1e-15)) const; + + // Return true if two polygons have boundaries that are within "max_error" + // of each other along their entire lengths. More precisely, there must be + // a bijection between the two sets of loops such that for each pair of + // loops, "a_loop->BoundaryNear(b_loop)" is true. + bool BoundaryNear(const S2Polygon& b, + S1Angle max_error = S1Angle::Radians(1e-15)) const; + + //////////////////////////////////////////////////////////////////////// + // S2Region interface (see s2region.h for details): + + // GetRectBound() returns essentially tight results, while GetCapBound() + // might have a lot of extra padding. Both bounds are conservative in that + // if the loop contains a point P, then the bound contains P also. + S2Polygon* Clone() const override; + S2Cap GetCapBound() const override; + S2LatLngRect GetRectBound() const override; + + bool Contains(const S2Cell& cell) const override; + bool MayIntersect(const S2Cell& cell) const override; + // Return true if the polygon contains the given point. Point containment + // is defined such that if the sphere is subdivided into polygons, every + // point is contained by exactly one polygons. This implies that polygons + // do not necessarily contain their vertices. + bool Contains(const S2Point& p) const override; + + // Encode the polygon with about 4 bytes per vertex, assuming the vertices + // have all been snapped to the centers of S2Cells at a given level + // (typically with InitToSnapped). The other vertices are stored using 24 + // bytes. Decoding a polygon encoded this way always returns the original + // polygon, without any loss of precision. + void Encode(Encoder* const encoder) const override; + + // Decode a polygon encoded with Encode(). + bool Decode(Decoder* const decoder) override; +}; +``` diff --git a/docs/devguide/coding_style.md b/docs/devguide/coding_style.md new file mode 100644 index 00000000..f059c767 --- /dev/null +++ b/docs/devguide/coding_style.md @@ -0,0 +1,923 @@ +--- +title: S2 C++ Style Guide +--- + + + +## Purpose + +This guide documents the C++ coding conventions used by the S2 Geometry Library. +It does not attempt to justify why these conventions are used, and indeed in +some cases they could be improved upon. + +In general, the S2 C++ coding style is built on the [Google C++ Style +Guide](https://google.github.io/styleguide/cppguide.html). However it's +important to remember that the S2 library has a long history, and many parts of +the library were developed before modern C++ features even existed. The S2 +library has many clients both within and outside of Google (through the open +source release), which can make it very difficult to adopt new conventions. + +**The overriding objective when adding new code should be to maintain +consistency within the library.** In cases where the S2 conventions conflict +with the current Google style guide, the existing S2 conventions should take +precedence. Changes in coding style should only be made when it is practical to +adopt them across the entire library, including the important step of updating +all affected clients so that any deprecated classes or methods can be removed. + +The first step when adding new classes or methods should be to find similar +classes or methods that already exist in the library, and then follow their +conventions as closely as possible. + +## External Dependencies + +The S2 library is built on many platforms, and has an open source version with +many external clients. New dependencies should be avoided. In particular, note +that S2 does not currently depend on the following: + + - GMock (except for a few files that not currently part of the open-source + release). Instead please limit yourself to the more basic features of GUnit. + + - Protocol buffers. + + - Any Google maps or geo library. + + - Many features of Abseil, including `StatusOr` (see "Error Handling" below). + +## Filenames + +The S2 library has a flat directory structure. Generally each header file +contains either a single class (plus any associated helper classes) or a +collection of functions in a namespace. + +Filenames should be in lowercase with words separated by underscores. Note that +"s2", "s1", etc. are not considered to be words. Examples: + + - `s2builder.h` + - `s2cell_id.h` + - `encoded_s2shape_index.h` + +For historical reasons, `latlng` is treated as one word even though the +corresponding classes are camel-cased. For example, `S2LatLngRect` is defined +in `s2latlng_rect.h`. + +Files that define functions or classes within the `s2shapeutil` and +`s2builderutil` namespaces include that namespace as a prefix: + + - `s2shapeutil_conversion.h` + - `s2shapeutil_count_edges.h` + - `s2builderutil_lax_polyon_layer.h` + - `s2builderutil_snap_functions.h` + +Classes or functions that are needed in multiple places within the S2 library +but are not part of the public API can be put in header files with the +`_internal.h` suffix. Examples: + + - `s2coords_internal.h` + - `s2predicates_internal.h` + +This can help remove clutter from the public API and also facilitate testing. +Note that the internal content can be implemented in the same `.cc` file as the +public part of the API. + +"Util" header files are strongly discouraged. Instead try to break your API +into meaningful components. The only remaining "util"-style header file in the +S2 library is `s2pointutil.h` which contains about ten functions for +manipulating points. + +## Documentation + +### Comment Style + +Please put effort into your comments; they are the most important part of your +code. + +Comments should be full sentences, starting with a capital letter and ending +with a period. This rule applies even to partial-line comments when reasonable: + +```c++ + enum IndexStatus { + STALE, // There are pending updates. + UPDATING, // Updates are currently being applied. + FRESH, // There are no pending updates. + }; +``` + +For bonus points, put two spaces after your periods: contrary to popular belief, +this is not some archaic convention from the days of typewriters, but rather an +effort to mimic professional typography (take a close look at a book sometime). + +### Header Files + +Public classes and methods should be thoroughly documented in the header file. +Important topics to cover include: + + - Description of the problem being solved + - Code samples showing how to use the API + - Limitations and known problems + +Requirements of classes and methods may be documented via a "REQUIRES" comment. +For example: + +```c++ +// Encodes an unsigned integer in little-endian format using "length" bytes. +// (The client must ensure that the encoder's buffer is large enough.) +// +// REQUIRES: T is an unsigned integer type. +// REQUIRES: 2 <= sizeof(T) <= 8 +// REQUIRES: 0 <= length <= sizeof(T) +// REQUIRES: value < 256 ** length +// REQUIRES: encoder->avail() >= length +template +void EncodeUintWithLength(T value, int length, Encoder* encoder); +``` + +The following callouts may also be useful: + + - `ENSURES:` properties guaranteed by the class/method. + - `DEFAULT:` documents the default value of an option. + - `CAVEAT:` an unexpected property or potential pitfall. + - `NOTE:` to draw special attention to a comment. + +### Implementations + +Implementation details should be documented where the corresponding types or +methods are defined, not where they are declared. This means that types and +member variables declared in the header file should be documented there, whereas +private non-inline methods should be documented in the .cc file. This +convention makes it easier to find documentation when looking at code, and also +makes it more likely that the documentation will be updated when implementations +are changed. + +Here is an example of documenting a private method in a `.cc` file: + +```c++ +// Outputs the current buffered path (which is assumed to be a loop), and +// resets the state to prepare for buffering a new loop. +void S2BufferOperation::OutputPath() { + op_.AddLoop(path_); + path_.clear(); // Does not change capacity. + have_input_start_ = false; + have_offset_start_ = false; +} +``` + +### External Documentation + +In addition to header file comments, significant new classes should have +documentation added in the "g3doc" directory. Typically this will expand +significantly on what is in the header file, giving more background information +and examples. Documentation added here will be published on the `s2geometry.io` +website. + +Ideally documentation should be added and/or updated in the same changelist that +adds the relevant code. + +## Namespaces + +The S2 library as a whole is not defined within a namespace. (Namespaces did +not even exist when it was first designed.) Instead most of the global symbols +include `S2` in their names (except for a few classes that use `S1`, `R2`, +etc. instead.) + +Nevertheless, some parts of the library do use namespaces. In general the names +of namespaces should be in lowercase with words separated by underscores +(e.g. `s2polyline_alignment`). Note that the "util" suffix is not considered a +word, so we have: + + - `namespace s2shapeutil { ... }` + - `namespace s2builderutil { ... }` + +There is also a namespace called `S2` that contains many of the low-level +functions defined by the library (e.g. `S2::TurnAngle`, `S2::RobustCrossProd`, +etc.) Note that while class names often include `S2` as a prefix rather +than being placed in a namespace, this convention is not typically used for +global functions or constants and so these entities should be put in a namespace +instead. + +In general S2 does not use nested namespaces, except for internal content which +may be placed in a nested namespace called `internal`. This technique can be +used alone to protect internal content (e.g. the file `s2testing.h` defines some +content in the `S2::internal` namespace) and/or with an `_internal.h` header +file. + +## Using Declarations + +`using` declarations are used extensively in `.cc` files in order to reduce code +verbosity. This includes types and functions from the `std` and `absl` +namespaces as well as aliases for nested types. For example, +`s2buffer_operation_test.cc` has the following declarations: + +```c++ +using absl::make_unique; +using std::max; +using std::string; +using std::unique_ptr; +using std::vector; + +using EndCapStyle = S2BufferOperation::EndCapStyle; +using PolylineSide = S2BufferOperation::PolylineSide; +``` + +Note that this is allowable only in `.cc` files, not header files, and that +`using namespace` declarations are prohibited. + +## Classes + +Class names in S2 should be nouns with the first letter of each word +capitalized. Examples: + + - `S2Shape` + - `S2LatLngRect` + - `S1Angle` + +Short, descriptive names are highly preferred. It is reasonable to put as much +thought into the names of your classes, methods, variables, etc. as you put into +the code itself. Good names are essential to making APIs easier to understand. + +The names of subtypes should be obviously related to their parent class. For +example, `S2LaxPolygonShape` is a subtype of `S2Shape`, and +`s2builderutil::S2CellIdSnapFunction` is a subtype of `S2Builder::SnapFunction`. + +### Index, Query, Operation + +Classes in S2 should be designed so as to separate the data representation from +the actual operations as much as possible.A class that contains many geometric +objects for the purpose of performing operations on them should preferentially +be called an `Index`. Such classes should define only the minimum number of +functions necessary to access the data in the index. Examples: + + - `S2ShapeIndex` + - `S2CellIndex` + - `S2PointIndex` + +Rather than having a single class with a wide API, prefer to have several +classes or standalone functions with narrow APIs. For example, rather than +having a method in `S2ShapeIndex` to count the total number of edges, instead +there is a separate function to do this defined in `s2shapeutil_count_edges.h`. + +A class that does something useful with an index should preferentially be called +a `Query`, unless it can be viewed as function that maps an index to another +similar index in which case it may be called an `Operation`. Examples: + + - `S2ClosestPointQuery` + - `S2ConvexHullQuery` + - `S2CrossingEdgeQuery` + - `S2BooleanOperation` + - `S2BufferOperation` + +`Index` classes should be thread-safe but `Query` and `Operation` classes +generally are not. Instead the expectation is that each thread should create +its own `Query` objects (typically on the stack) as needed. The fact that +`Query` classes are not thread-safe allows them to maintain internal state while +executing the query and/or to cache information between queries. + +`Query` and `Operation` classes typically follow one of two design patterns. +One is the builder pattern, where data can be added incrementally using several +method calls, and then a different method is called to compute the result. +Examples following this pattern include `S2WindingOperation`, +`S2BufferOperation`, `S2ConvexHullQuery`, etc. + +The other common pattern is a pure function pattern where the arguments are +supplied and the result is returned during a single function call. Examples of +this pattern include `S2ClosestEdgeQuery`, `S2ClosestPointQuery`, +`S2ContainsPointQuery`, etc. Note that even though this pattern in theory +allows the query method to be stateless and therefore thread-safe, S2 classes +avoid this in favor of rather allowing query objects to maintain whatever state +they may need in order to make queries as efficient as possible. This means +that **query methods should not be const** even if the current implementation of +the method would allow this, since this might prevent future optimizations. + +When query methods need to return several values of different types, the query +method should return a result object. For example, `S2ClosestEdgeQuery` returns +a `Result` object with methods such as `distance()`, `shape_id()`, and +`edge_id()`. + +This technique is strongly preferred over the alternative technique of having +the query object store the latest result in its own state and having methods to +retrieve the result field-by-field. For example, a Result object makes it much +easier to combine the results of several queries: + +``` + auto a_result = query.GetResult(a); + auto b_result = query.GetResult(b); + ... do something with both results ... +``` + +Note that as long as the optimizer does its job correctly, returning a result by +value does not involve any actual copying. Instead the caller passes a hidden +pointer and the result is constructed directly on the caller's stack. + +### Options + +Classes with options should have a nested `Options` class. For +example, `S2Builder` has an `S2Builder::Options` class. Options +classes generally consist of a set of fields with "get" and "set" methods. For +example: + +```c++ +class S2RegionCoverer { + public: + class Options { + public: + Options(); + + int max_cells() const; + void set_max_cells(int max_cells); + ... + }; + explicit S2RegionCoverer(const Options& options); + const Options& options() const; + ... +}; +``` + +The set methods generally return `void`, but it is also acceptable to return a +reference to the Options object in order to allow chained initialization: + +``` + S2Foo foo{S2Foo::Options().set_bar(x).set_baz(y)}; +``` + +Options are usually read-only once the relevant object has been constructed. +Classes that allow options to be modified after construction do so via a +`mutable_options()` method: + +```c++ + Options* mutable_options(); +``` + +Options classes are copyable and movable. + +## Methods + +Accessor method names are nouns, written in lowercase with words separated by +underscores. Examples: + + - `S1Angle::degrees()` + - `S2Builder::options()` + +Non-accessor method names should be verbs, written in uppercase with the first +letter of each word capitalized. Examples: + + - `S2Builder::StartLayer(...)` + - `S2Builder::Build(...)` + - `S2ShapeIndex::Minimize()` + +Methods that return something but that don't have an obvious verb should +typically start with `Get`. For example: + + - `S2ShapeIndex::GetShape(int id)` + - `S2::GetPointOnLine(...)` + +If the method may return an empty set of results, the verb `Find` can also be +used. For example: + + - `S2ClosestEdgeQuery::FindClosestEdges(...)` + +There are exceptions for well-known geometric operations that yield objects of +the same type. For example, `S2Cap::Union` computes the union of two `S2Cap` +objects. This naming convention is generally reserved for low-level geometric +objects such as `S1Interval`, `S2LatLngRect`, `S2Cap`, etc. + +Methods that return boolean values (predicates) should be active verbs +(e.g. `S2LatLngRect::Contains()`) or start with the word "is" +(e.g. `S2Shape::is_empty()`, `S2::IsUnitLength()`). (Note that this convention +differs from `std::vector::empty()` for example.) + +### Inline Methods + +Except for accessor methods, inline methods should be declared and defined +separately (just like non-inline methods). Note that only the definition needs +to be declared inline, not the declaration. Inline definitions should be placed +at the bottom of the corresponding header file after a comment saying +"Implementation details follow". For example: + +```c++ +class S2PaddedCell { + public: + // Construct an S2PaddedCell for the given cell id and padding. + S2PaddedCell(S2CellId id, double padding); + + S2CellId id() const { return id_; } + double padding() const { return padding_; } + int level() const { return level_; } + + // Return the bound for this cell (including padding). + const R2Rect& bound() const { return bound_; } + + // Return the (i,j) coordinates for the child cell at the given traversal + // position. The traversal position corresponds to the order in which child + // cells are visited by the Hilbert curve. + void GetChildIJ(int pos, int* i, int* j) const; + + ... +}; + + +////////////////// Implementation details follow //////////////////// + + +inline void S2PaddedCell::GetChildIJ(int pos, int* i, int* j) const { + int ij = S2::internal::kPosToIJ[orientation_][pos]; + *i = ij >> 1; + *j = ij & 1; +} +``` + +### Constructors, Init, Factory Methods + +Virtually all classes should have constructors. If the constructor parameters +may be ambiguous, static factory methods should be used instead. For example, +`S2LatLng` has the following factory methods: + +```c++ + static constexpr S2LatLng FromRadians(double lat_radians, double lng_radians); + static constexpr S2LatLng FromDegrees(double lat_degrees, double lng_degrees); + static constexpr S2LatLng FromE5(int32_t lat_e5, int32_t lng_e5); + static constexpr S2LatLng FromE6(int32_t lat_e6, int32_t lng_e6); + static constexpr S2LatLng FromE7(int32_t lat_e7, int32_t lnt_e7); +``` + +Factory method names should begin with "From", except for simple and frequently +used classes where brevity takes precedence. For example, `S1Angle` has the +factory methods `S1Angle::Radians(double radians)` and `S1Angle::Degrees(double +degrees)`. + +Factory methods should also be used when the parameter order may be ambiguous. +Examples: + +```c++ + static R2Rect FromCenterSize(const R2Point& center, const R2Point& size); + static S2CellId FromFacePosLevel(int face, uint64 pos, int level); +``` + +Complex classes that allocate memory should also have a default constructor and +an `Init` method. This allows clients to declare instances of these classes as +member variables and defer their initialization as necessary. This allows +clients to avoid allocating objects on the heap and accessing them through +`unique_ptr<>` just because the arguments needed to initialize them are not +readily available in their own constructor. This not only reduces memory +allocation but also increases locality of reference. + +So for example, `S2Builder` has the following: + +```c++ +class S2Builder { + ... + // Default constructor; requires Init() to be called. + S2Builder(); + + // Convenience constructor that calls Init(). + explicit S2Builder(const Options& options); + + // Initializes an S2Builder with the given options. + void Init(const Options& options); + + const Options& options() const { return options_; } + ... +}; +``` + +Factory methods that return `unique_ptr` should be avoided, since this +forces clients to allocate memory and also reduces locality of reference. + +## Types + +### Templates + +Templated types and functions are discouraged in the public parts of the API. +Templated code is necessarily harder to document and understand, and the +generalization that it offers is frequently overrated and underused. + +Templates may be used internally when this substantially reduces code +duplication or otherwise makes the library easier to maintain. For example, the +S2 classes `S2ClosestEdgeQuery` and `S2FurtherEdgeQuery` are each wrappers +around a templated internal class called `S2ClosestEdgeQueryBase`. + +### Distances and Angles + +Distances and angles should be represented using the classses `S1Angle` or +`S1ChordAngle`. S2 APIs should always use these classes in preference to +accepting distances in radians, degrees, meters, etc. since it avoids the need +to specify units. (Note that the S2 library operates on the unit sphere, not +the Earth's surface, and therefore the natural unit of distance is the radian.) + +Whenever possible, classes should support both `S1Angle` and `S1ChordAngle`. +`S1ChordAngle` is frequently used internally, since (1) it is much faster to +compute and (2) it supports the exact distance predicates that are needed to +implement algorithms robustly (see `s2predicates.h`). Methods that accept an +`S1Angle` often convert it to `S1ChordAngle` internally. + +Note that `S1ChordAngle` is only capable of representing distances up to 180 +degrees, and that its accuracy declines for distances approaching 180 degrees. +Distances which may be greater than 180 degrees, such as polyline lengths or +polygon perimeters, should always be represented using `S1Angle`. + +### Parameters + +Small types (16 bytes or less) should generally be passed by value. This +includes the following commonly used S2 types: + + - `S2CellId` + - `S2LatLng` + - `S1Angle` + - `S1ChordAngle` + +Types larger than 16 bytes should be passed as follows: + + - If ownership is being transferred, use `unique_ptr`. + + - If the callee keeps a pointer to the object that outlives the function call, + use `const T*`. For example, `S2ClosestEdgeQuery` keeps a pointer to its + `S2ShapeIndex` constructor argument, so it accepts that argument by const + pointer: + + ```c++ + explicit S2ClosestEdgeQuery(const S2ShapeIndex* index, + const Options& options = Options()); + ``` + + - Otherwise, the parameter should be passed by const reference. Note in + particular that `S2Point` is passed by const reference. Const pointers + should not be used unless the callee keeps a pointer as outlined above. + + - The only exception to the rule above is when the callee needs to make a copy + of the argument **and*** the given type has an efficient move operator (i.e., + one that transfers ownership of data rather than simply copying it). In this + case the parameter may be passed by value rather than by const reference. + For example: + + ```c++ + class S2CellUnion final : public S2Region { + ... + + // Constructs a cell union with the given S2CellIds, then calls Normalize() + // to sort them, remove duplicates, and merge cells when possible. + // + // The argument is passed by value, so if you are passing a named variable + // and have no further use for it, consider using std::move(). + // + // A cell union containing a single S2CellId may be constructed like this: + // + // S2CellUnion example({cell_id}); + explicit S2CellUnion(std::vector cell_ids); + ... + } + + inline S2CellUnion::S2CellUnion(std::vector cell_ids) + : cell_ids_(std::move(cell_ids)) { + Normalize(); + } + ``` + + Note that this guidance is considerably more strict than advocated by the + current Google style guide. Do not pass parameters larger than 16 bytes + unless (1) the type has an efficient `std::move` operator and (2) the callee + makes a copy of the argument using `std::move`. + +#### Const Value Parameters + +Parameters passed by value should never be "const" in method declarations. They +should only be declared const (if desired) in method **definitions**. For +example: + +```c++ +class S2CellId { + ... + // Return true if the given cell is contained within this one. + bool contains(S2CellId other) const; + ... +}; + + +////////////////// Implementation details follow //////////////////// + + +inline bool S2CellId::contains(const S2CellId other) const { + DCHECK(is_valid()); + DCHECK(other.is_valid()); + return other >= range_min() && other <= range_max(); +} +``` + +Note that most S2 code does not bother declaring local variables and parameters +passed by value as "const" even when those values are not modified, since the +additional documentation benefit often does not outweigh the additional code +clutter. + +### Geometry Output + +Classes that construct geometry and return it to the client should do so via an +`S2Builder::Layer` parameter. For example, the `S2BufferOperation` constructor +looks like this: + +```c++ + explicit S2BufferOperation(std::unique_ptr result_layer, + const Options& options = Options{}); +``` + +When the client eventually calls the `S2BufferOperation::Build()` method, the +output geometry is sent to the given `result_layer`. This allows the client to +have complete control over how the output geometry is represented. There are +various possible layer types that build different representations such as +`S2Polyon` or `S2LaxPolygonShape`. + +### Obsolete Types + +Clients are discouraged from using certain S2 classes for new code. Classes +that have officially been deprecated are marked as such using the +`ABSL_DEPRECATED` macro. Currently the only such class is `S2PolygonBuilder` +which has been superseded by the far more robust and flexible `S2Builder`. + +However it is important to note that the classes `S2Polygon` and `S2Polyline` +are also obsolete and should ideally be deprecated at some point in the future. +New code should prefer to use their `S2Shape`-based replacements, +`S2LaxPolygonShape` and `S2LaxPolylineShape`. These new classes can be +constructed and decoded far more efficiently than `S2Polygon` and `S2Polyline`, +and are also capable of representing degenerate geometry if desired (which is +what the `Lax` in their names refers to). + +## Error Handling + +In general S2 code should try to avoid generating runtime errors. + +### Parameter Range Checking + +If a parameter has a limited valid range, the preferred technique in S2 is to +`DCHECK` that the parameter is valid and then clamp it to the +allowed range. For example: + +```c++ + class S2BufferOperation::Options { + ... + // Specifies the allowable error when buffering, expressed as a fraction + // of buffer_radius(). + // + // REQUIRES: error_fraction() >= kMinErrorFraction + // REQUIRES: error_fraction() <= 1.0 + // + // DEFAULT: 0.01 (i.e., maximum error of 1%) + static constexpr double kMinErrorFraction = 1e-6; + double error_fraction() const; + void set_error_fraction(double error_fraction); + ... + }; + +void S2BufferOperation::Options::set_error_fraction(double error_fraction) { + DCHECK_GE(error_fraction, kMinErrorFraction); + DCHECK_LE(error_fraction, 1.0); + error_fraction_ = max(kMinErrorFraction, min(1.0, error_fraction)); +} +``` + +### Invariants + +Internal invariants should be enforced via `DCHECK`. This is an important +aspect of both documentation and testing. For example, the following is a +private method that should only be called when certain conditions are true: + +```c++ +void S2BufferOperation::BufferEdgeAndVertex(const S2Point& a, const S2Point& b, + const S2Point& c) { + DCHECK_NE(a, b); + DCHECK_NE(b, c); + DCHECK_NE(buffer_sign_, 0); + if (!tracker_.ok()) return; + ... +} +``` + +Programming errors within the S2 library should not be converted to runtime +errors if at all possible. + + +### S2Error + +For historical reasons, S2 defines its own error class `S2Error` consisting of +an error code and a status message. This class is used exclusively within the +library in preference to other error representations such as `absl::Status`. +While it may be desirable to eventually convert the library to use a different +error representation, the overriding principle is to maintain consistency within +the library so that users (including external users of the open source version) +do not need to deal with two different error representations at once. + +Here is a partial definition of the `S2Error` class: + +```c++ +class S2Error { + public: + enum Code { + OK = 0, // No error. + + //////////////////////////////////////////////////////////////////// + // Generic errors, not specific to geometric objects: + + UNKNOWN = 1000, // Unknown error. + UNIMPLEMENTED = 1001, // Operation is not implemented. + OUT_OF_RANGE = 1002, // Argument is out of range. + INVALID_ARGUMENT = 1003, // Invalid argument (other than a range error). + FAILED_PRECONDITION = 1004, // Object is not in the required state. + INTERNAL = 1005, // An internal invariant has failed. + DATA_LOSS = 1006, // Data loss or corruption. + RESOURCE_EXHAUSTED = 1007, // A resource has been exhausted. + CANCELLED = 1008, // Operation was cancelled. + + ... + }; + S2Error() : code_(OK), text_() {} + + // Convenience constructor that calls Init(). + template + S2Error(Code code, const absl::FormatSpec& format, + const Args&... args); + } + + // Set the error to the given code and printf-style message. Note that you + // can prepend text to an existing error by calling Init() more than once: + // + // error->Init(error->code(), "Loop %d: %s", j, error->text().c_str()); + template + void Init(Code code, const absl::FormatSpec& format, + const Args&... args); + + bool ok() const { return code_ == OK; } + Code code() const { return code_; } + std::string text() const { return text_; } + + // Clear the error to contain the OK code and no error message. + void Clear(); + + private: + Code code_; + std::string text_; +}; +``` + +Errors are typically initialized using `S2Error::Init()`, although they can also +be constructed and copied. For example: + +```c++ + error->Init(S2Error::RESOURCE_EXHAUSTED, + "Memory limit exceeded (tracked usage %d bytes, limit %d bytes)", + usage_bytes_, limit_bytes_); +``` + +`S2Error` is typically passed as a pointer parameter and is required to be +non-`nullptr`. (Note that the current Google style recommendation is to pass such +parameters by reference, which precludes `nullptr` arguments at the expense of +giving up the visible documentation hint that the argument may be modified.) +Public methods typically call `error->Clear()` to remove any existing error, +whereas private methods either set the error or leave its state unchanged. + +Methods that would otherwise return `void` often return `error->ok()` instead in +order to simplify client code. For example: + +```c++ +class S2BooleanOperation { + ... + // Executes the given operation. Returns true on success, and otherwise + // sets "error" appropriately. (This class does not generate any errors + // itself, but the S2Builder::Layer might.) + bool Build(const S2ShapeIndex& a, const S2ShapeIndex& b, S2Error* error); +} +``` + +Client code might then look like this: + +```c++ + bool DoSomething(..., S2Error* error) { + ... + if (!op.Build(a, b, error)) return false; + ... + } +``` + +If a function has a void or other non-boolean result, the caller must check +`error->ok()` explicitly: + +```c++ + void DoSomething(S2Builder::Layer* layer, const S2Builder::Graph& graph, + S2Error* error) { + layer->Build(graph, error); + if (!error->ok()) return; + ... + } +``` + +## Testing + +Tests should attempt to achieve excellent coverage with a minimum of code. This +will usually mean writing helper functions and/or test fixtures in order to +avoid code duplication. Do not be afraid to write test code with loops, +multiple levels of helpers, etc., if this will result in better testing with +less code duplication. + +S2 tests should always be run **in a debug build** (not the default "fast +build") since this will give full line number information on all your stack +traces. In other words, build S2 with `-c dbg` when you are testing new code. + +Non-optimized builds also enable various internal consistency checks. This +includes not only the `DCHECK` macros but also the `--s2debug` flag which +enables certain internal validity checking. (This flag can also be enabled by +hand if you are suspect a problem with geometry that is too large or slow to +process using a debug build.) + +Tests should focus on the functionality provided by your class, not the +functionality of the classes that your class depends on. For example, if your +class uses `S2Builder` then you should not write tests whose purpose is simply +to verify the guarantees that `S2Builder` makes (which are thoroughly tested +elsewhere). + +### Allow For Errors + +Tests should generally allow for small amounts of error in their results rather +than expecting an exact match. This helps to make tests less fragile, avoiding +the need to update them whenever small changes to the underlying algorithms are +made. + +Useful Gunit helper functions include `EXPECT_DOUBLE_EQ` and `EXPECT_NEAR`. +It can also be useful to snap results to a fixed resolution using +`s2builderutil::IntLatLngSnapFunction()`. + +"Golden" files (which compare the results of a series of operations to a +"golden" result deemed to be correct) are strictly prohibited.` Such tests are +inherently fragile and will break whenever the smallest details of S2 algorithms +are changed. + +### S2Testing + +The header file `s2testing.h` provides many utility functions and classes that +are useful when writing tests. In addition the distance-related functions +mentioned above, there is extensive support for generating random geometry of +various kinds. + +Support code for writing tests can also be found in the following files: + + - `thread_testing.h` + - `s2shapeutil_testing.h` + - `s2closest_edge_query_testing.h` + +### S2Earth + +S2 tests are not allowed to depend on `s2earth.h` (which provides basic +functionality for converting distances on the Earth's surface to S2-compatible +types). Instead tests should use the even more basic functionality provided in +`s2testing.h`. For example: + +```c++ + S1Angle tolerance = S2Testing::MetersToAngle(1); +``` + +### Randomization and Robustness + +S2 makes heavy use of randomized tests. This is particularly useful as a way to +test robustness guarantees. Typically, a specialized function is written to +generate test cases that are very likely to trigger problems, and then a large +number of cases are checked for correctness. Here is an example: + +```c++ +// Chooses a random S2Point that is often near the intersection of one of the +// coodinates planes or coordinate axes with the unit sphere. (It is possible +// to represent very small perturbations near such points.) +S2Point ChoosePoint() { + S2Point x = S2Testing::RandomPoint(); + for (int i = 0; i < 3; ++i) { + if (S2Testing::rnd.OneIn(3)) { + x[i] *= pow(1e-50, S2Testing::rnd.RandDouble()); + } + } + return x.Normalize(); +} + +TEST(S2, ProjectError) { + for (int iter = 0; iter < 1000; ++iter) { + SCOPED_TRACE(StrCat("Iteration ", iter)); + S2Testing::rnd.Reset(iter); // Easier to reproduce a specific case. + S2Point a = ChoosePoint(); + S2Point b = ChoosePoint(); + S2Point n = S2::RobustCrossProd(a, b).Normalize(); + S2Point x = S2Testing::SamplePoint(S2Cap(n, S1Angle::Radians(1e-15))); + S2Point p = S2::Project(x, a, b); + EXPECT_LT(s2pred::CompareEdgeDistance( + p, a, b, S1ChordAngle(S2::kProjectPerpendicularError)), 0); + } +} +``` + +Note the use of `SCOPED_TRACE` so that any failing iterations can be identified. +This is obviously only useful if the random number seed is deterministic, which +is one reason that `S2Testing` defines its own random number generator available +as `S2Testing::rnd`. Tests **should not** use random number generators with +non-deterministic seeds, since this makes problems impossible to debug. + +Also note the use of `S2Testing::rnd.Reset(iter)` which allows a specific +failing test case to be reproduced more easily by simply replacing 'iter' by the +failing iteration number. (Note that most S2 code actually calls +`rnd.Reset(iter + 1)` because it turns out that `rnd.Reset(0)` and +`rnd.Reset(1)` have the same effect.) + +Observe the use of `S2Testing` helper methods such as `RandomPoint` and +`SamplePoint`. Also note the use of `pow(10**x, S2Testing::rnd.RandDouble())`, +which produces a value whose *exponent* is uniformly distributed between "x" +and 0. This is a useful technique for generating numbers that are +well-distributed over a wide range of magnitudes. diff --git a/docs/devguide/cpp/examples/CMakeLists.txt b/docs/devguide/cpp/examples/CMakeLists.txt new file mode 100644 index 00000000..d5858e01 --- /dev/null +++ b/docs/devguide/cpp/examples/CMakeLists.txt @@ -0,0 +1,4 @@ +add_executable(point_index point_index.cc) +target_link_libraries(point_index LINK_PUBLIC s2testing s2) +add_executable(term_index term_index.cc) +target_link_libraries(term_index LINK_PUBLIC s2testing s2) diff --git a/docs/devguide/cpp/examples/point_index.cc b/docs/devguide/cpp/examples/point_index.cc new file mode 100644 index 00000000..b924055c --- /dev/null +++ b/docs/devguide/cpp/examples/point_index.cc @@ -0,0 +1,56 @@ +// Copyright 2017 Google Inc. All Rights Reserved. +// Author: ericv@google.com (Eric Veach) +// +// This example shows how to build and query an in-memory index of points +// using S2PointIndex. + +#include +#include +#include + +#include "base/commandlineflags.h" +#include "base/init_google.h" // MOE:strip_line +#include "geostore/geomodel/public/s2earth.h" +#include "third_party/absl/flags/flag.h" +#include "util/geometry/s1chord_angle.h" +#include "util/geometry/s2closest_point_query.h" +#include "util/geometry/s2point_index.h" +#include "util/geometry/s2testing.h" + +DEFINE_int32(num_index_points, 10000, "Number of points to index"); +DEFINE_int32(num_queries, 10000, "Number of queries"); +DEFINE_double(query_radius_km, 100, "Query radius in kilometers"); + +// MOE:begin_strip +using geostore::S2Earth; + +// MOE:end_strip +int main(int argc, char **argv) { + // MOE:begin_strip + // TODO(jrosenstock,b/209396162): We should call + // gflags::ParseCommandLineFlags here, but that would break WITH_FLAGS=off. + InitGoogle(argv[0], &argc, &argv, true); + // MOE:end_strip + // Build an index containing random points anywhere on the Earth. + S2PointIndex index; + for (int i = 0; i < absl::GetFlag(FLAGS_num_index_points); ++i) { + index.Add(S2Testing::RandomPoint(), i); + } + + // Create a query to search within the given radius of a target point. + S2ClosestPointQuery query(&index); + query.mutable_options()->set_max_distance(S1Angle::Radians( + S2Earth::KmToRadians(absl::GetFlag(FLAGS_query_radius_km)))); + + // Repeatedly choose a random target point, and count how many index points + // are within the given radius of that point. + int64_t num_found = 0; + for (int i = 0; i < absl::GetFlag(FLAGS_num_queries); ++i) { + S2ClosestPointQuery::PointTarget target(S2Testing::RandomPoint()); + num_found += query.FindClosestPoints(&target).size(); + } + + std::printf("Found %" PRId64 " points in %d queries\n", num_found, + absl::GetFlag(FLAGS_num_queries)); + return 0; +} diff --git a/docs/devguide/cpp/examples/term_index.cc b/docs/devguide/cpp/examples/term_index.cc new file mode 100644 index 00000000..d71a1226 --- /dev/null +++ b/docs/devguide/cpp/examples/term_index.cc @@ -0,0 +1,112 @@ +// Copyright 2017 Google Inc. All Rights Reserved. +// Author: ericv@google.com (Eric Veach) +// +// This example shows how to add spatial data to an information retrieval +// system. Such systems work by converting documents into a collection of +// "index terms" (e.g., representing words or phrases), and then building an +// "inverted index" that maps each term to a list of documents (and document +// positions) where that term occurs. +// +// This example shows how to convert spatial data into index terms, which can +// then be indexed along with the other document information. + +#include +#include +#include +#include +#include +#include + +#include "base/commandlineflags.h" +#include "base/init_google.h" // MOE:strip_line +#include "geostore/geomodel/public/s2earth.h" +#include "third_party/absl/container/btree_set.h" +#include "third_party/absl/container/flat_hash_map.h" +#include "third_party/absl/flags/flag.h" +#include "util/geometry/s2cap.h" +#include "util/geometry/s2point_index.h" +#include "util/geometry/s2region_term_indexer.h" +#include "util/geometry/s2testing.h" + +DEFINE_int32(num_documents, 10000, "Number of documents"); +DEFINE_int32(num_queries, 10000, "Number of queries"); +DEFINE_double(query_radius_km, 100, "Query radius in kilometers"); + +// MOE:begin_strip +using geostore::S2Earth; + +// MOE:end_strip +// A prefix added to spatial terms to distinguish them from other index terms +// (e.g. representing words or phrases). +static const char kPrefix[] = "s2:"; + +int main(int argc, char** argv) { + // MOE:begin_strip + // TODO(jrosenstock,b/209396162): We should call + // gflags::ParseCommandLineFlags here, but that would break WITH_FLAGS=off. + InitGoogle(argv[0], &argc, &argv, true); + // MOE:end_strip + // Create a set of "documents" to be indexed. Each document consists of a + // single point. (You can easily substitute any S2Region type here, or even + // index a mixture of region types using std::unique_ptr. Other + // region types include polygons, polylines, rectangles, discs, buffered + // geometry, etc.) + std::vector documents; + documents.reserve(absl::GetFlag(FLAGS_num_documents)); + for (int docid = 0; docid < absl::GetFlag(FLAGS_num_documents); ++docid) { + documents.push_back(S2Testing::RandomPoint()); + } + + // We use a hash map as our inverted index. The key is an index term, and + // the value is the set of "document ids" where this index term is present. + absl::flat_hash_map> index; + + // Create an indexer suitable for an index that contains points only. + // (You may also want to adjust min_level() or max_level() if you plan + // on querying very large or very small regions.) + S2RegionTermIndexer::Options options; + options.set_index_contains_points_only(true); + S2RegionTermIndexer indexer(options); + + // Add the documents to the index. + for (int docid = 0; docid < documents.size(); ++docid) { + S2Point index_region = documents[docid]; + for (const auto& term : indexer.GetIndexTerms(index_region, kPrefix)) { + index[term].push_back(docid); + } + } + + // Convert the query radius to an angle representation. + S1Angle radius = S1Angle::Radians( + S2Earth::KmToRadians(absl::GetFlag(FLAGS_query_radius_km))); + + // Count the number of documents (points) found in all queries. + int64_t num_found = 0; + for (int i = 0; i < absl::GetFlag(FLAGS_num_queries); ++i) { + // Choose a random center for query. + S2Cap query_region(S2Testing::RandomPoint(), radius); + + // Convert the query region to a set of terms, and compute the union of + // the document ids associated with those terms. (An actual information + // retrieval system would do something more sophisticated.) + absl::btree_set candidates; + for (const auto& term : indexer.GetQueryTerms(query_region, kPrefix)) { + candidates.insert(index[term].begin(), index[term].end()); + } + + // "candidates" now contains all documents that intersect the query + // region, along with some documents that nearly intersect it. We can + // prune the results by retrieving the original "document" and checking + // the distance more precisely. + std::vector result; + for (int docid : candidates) { + if (!query_region.Contains(documents[docid])) continue; + result.push_back(docid); + } + // Now do something with the results (in this example we just count them). + num_found += result.size(); + } + std::printf("Found %" PRId64 " points in %d queries\n", num_found, + absl::GetFlag(FLAGS_num_queries)); + return 0; +} diff --git a/docs/devguide/cpp/quickstart.md b/docs/devguide/cpp/quickstart.md new file mode 100644 index 00000000..1bea5896 --- /dev/null +++ b/docs/devguide/cpp/quickstart.md @@ -0,0 +1,360 @@ +--- +title: Quick Start +--- + +Suppose that we have a few terabytes (or petabytes) of geographic data, and we +want to be able to query it efficiently. This brief tutorial shows how the S2 +library is useful for solving this problem. + +## Prerequisites + +Before going through this Quickstart, first makes sure your development +environment satisfies our +[platform requirements](/about/platforms) and +follow the installation instructions in the +[S2 Install Guide](/about/install). + +## S2PointIndex + +To keep things simple, let's suppose that our input data consists only of points +(rather than polylines, polygons, etc), and furthermore let's suppose that we +have only a few million points, so that they all fit in memory. (We will look +at the more general case below.) + +The easiest way to build an in-memory point index is to use an `S2PointIndex`, +as follows: + + // Build an index containing random points anywhere on the Earth. + S2PointIndex index; + for (int i = 0; i < FLAGS_num_index_points; ++i) { + S2Point point = S2Testing::RandomPoint(); + index.Add(point, i); + } + +This example uses `S2Testing::RandomPoint()` to generate points that are +randomly distributed over the Earth. Note that in the S2 library, [points are +represented as unit-length +vectors](/about/overview#unit-vectors-vs-latitudelongitude-pairs) +(the [`S2Point`](/devguide/basic_types#s2point) class), +rather than as latitude-longitude pairs. You can convert latitude-longitude +pairs to `S2Point` values using the +[`S2LatLng`](/devguide/basic_types#s2latlng) class, like +this: + + S2Point point(S2LatLng::FromDegrees(lat_degrees, lng_degrees)); + +Also notice that each point in the index has some auxiliary data associated +with it (in this case, the integer `i`). You can use this feature to map the +points back to other data structures. + +## S2ClosestPointQuery + +Now that we have our index, how do we query it? For example, suppose that we +want to find all the points within 100km of a target point. This can be done +using `S2ClosestPointQuery`, which provides methods that find the closest +point(s) to a given target object (such as a point, polyline, or polygon). +You can either find the *k* closest points to the target, or all points within +a given radius of the target, or both. + +In this example, we will find all points within 100km of a point target. +First we construct the `S2ClosestPointQuery` object: + + S2ClosestPointQuery query(&index); + query.mutable_options()->set_max_distance( + S1Angle::Radians(S2Earth::KmToRadians(FLAGS_query_radius_km))); + +Notice that we need to convert the radius in kilometers to an `S1Angle`. This +is because S2 models the Earth as *unit sphere* (a sphere with radius 1), +which means that the distance between two points is the same as the angle +between those points (in radians) measured from the sphere's center. S2 +provides helper function for converting distances to angles, such as the +`S2Earth::KmToRadians` function used above. (Note that because the Earth is +not quite spherical, distances computed in this way can have errors of up to +0.56%. You can use a geodesy library to compute distances more accurately, +but spherical distances are perfectly adequate for many purposes, such as +finding nearby restaurants.) + +Finally, we can find the points within the given radius of our target point +like this: + + S2ClosestPointQuery::PointTarget target(S2Testing::RandomPoint()); + auto results = query.FindClosestPoints(&target); + +This returns a vector of `Result` objects, where each result includes the +following fields: + + S2Point point() const; // The indexed point. + Data data() const; // The auxiliary data for this point. + S1ChordAngle distance() const; // The distance to this point. + +In our case `data()` is the integer value that we passed to `index.Add()` +above. + +Putting this all together, our program for indexing and querying points looks +like this (see `examples/point_index.cc`): + +```c++ + // Build an index containing random points anywhere on the Earth. + S2PointIndex index; + for (int i = 0; i < FLAGS_num_index_points; ++i) { + index.Add(S2Testing::RandomPoint(), i); + } + + // Create a query to search within the given radius of a target point. + S2ClosestPointQuery query(&index); + query.mutable_options()->set_max_distance( + S1Angle::Radians(S2Earth::KmToRadians(FLAGS_query_radius_km))); + + // Repeatedly choose a random target point, and count how many index points + // are within the given radius of that point. + int64 num_found = 0; + for (int i = 0; i < FLAGS_num_queries; ++i) { + S2ClosestPointQuery::PointTarget target(S2Testing::RandomPoint()); + num_found += query.FindClosestPoints(&target).size(); + } + + printf("Found %lld points in %d queries\n", num_found, FLAGS_num_queries); +``` + +## Indexing for Search + +Let's suppose again that we have a lot of geographic data, and we want to +query it efficiently along with other attributes (such as names and +descriptions) using an *information retrieval system*. Such systems index a +collection of *documents*, where each document consists of text and other +attributes. Each document is converted into a collection of *index terms* +(e.g., significant words or phrases), which are gathered into an *inverted +index* that maps each term to a list of documents where that term occurs. +(For a comprehensive introduction to information retrieval, see [Introduction +to Information Retrieval](https://nlp.stanford.edu/IR-book/){:target="_blank"} +by Christopher D. Manning et al., Cambridge University Press, 2008.) + +Our goal is not to solve the information retrieval problem in this tutorial, +but rather just to show how to add spatial data to an existing system. S2 +defines a special class called `S2RegionTermIndexer` that is designed to deal +with the problem of converting spatial data into index terms, which can then +be indexed along with the other document information. (In fact the method +described here can be used with any key-value storage system.) + +To keep things simple, let's assume that each document consists of a single +point (with no other information). Each document is assigned a "document id" +that can be used to retrieve it later: + + std::vector documents; + for (int docid = 0; docid < FLAGS_num_documents; ++docid) { + documents.push_back(S2Testing::RandomPoint()); + } + +As a substitute for the "inverted index" of an information retrieval system, +we will use a hash map. The key of the hash map is an index term, and the +value is the set of document ids where this term is present: + + std::unordered_map> index; + +We can now use `S2RegionTermIndexer` to index these documents as follows: + + // Create an indexer with appropriate options. + S2RegionTermIndexer::Options options; + S2RegionTermIndexer indexer(options); + + // Add the documents to the index. + static const char kPrefix[] = "s2:"; + for (int docid = 0; docid < documents.size(); ++docid) { + S2Point index_region = documents[docid]; + for (const auto& term : indexer.GetIndexTerms(index_region, kPrefix)) { + index[term].push_back(docid); + } + } + +For each document, the indexer returns a set of terms that should be +associated with the document for later retrieval. Each term is simply an +alphanumeric string. To distinguish these spatial indexing terms from other +index terms (such as words), the interface allows a prefix to be added (in +this case "s2:"). For example, an index term might look like "s2:34f14a9". + +## Querying for Search + +Let's now examine how to query our inverted index. Again, for simplicity +let's assume that we want to retrieve all documents that are within a given +radius of a test point. As before, the first step is to convert the query +radius to an `S1Angle`: + + S1Angle radius = + S1Angle::Radians(S2Earth::KmToRadians(FLAGS_query_radius_km)); + +Our query region is a disc centered around a target point. On the sphere, a +disc-shaped region such as this one is called a *spherical cap* +([`S2Cap`](/devguide/basic_types#s2cap)): + + S2Cap query_region(S2Testing::RandomPoint(), radius); + +Now we convert the query region to a set of *query terms*. The query terms +are constructed such that if the query region intersects the document region, +then the query terms and index terms are guaranteed to have at least one value +in common. This means that to find the set of intersecting documents, we +simply need to look up the documents associated with each query term and +compute their union. (Actual information retrieval systems do something more +sophisticated, but that doesn't concern us here.) + + std::set candidates; + for (const auto& term : indexer.GetQueryTerms(query_region, kPrefix)) { + candidates.insert(index[term].begin(), index[term].end()); + } + +"candidates" now contains all documents that intersect the query region, along +with some documents that nearly intersect it. At this point, an information +retrieval system would "score" documents by examining their content more +closely. In our case, for example, we might want to filter the candidates by +retrieving the original "document" and checking the distance to the target +more precisely. We can do this as follows: + + std::vector result; + for (int docid : candidates) { + if (!query_region.Contains(documents[docid])) continue; + result.push_back(docid); + } + +Now the `results` vector contains exactly the set of input documents that +intersect the query region. + +## Optimization for Points + +One final detail: when a spatial index contains only points, like this one, it +turns out that fewer query terms are needed than when the index also contains +non-point regions (such as polylines or polygons). You can take advantage of +this optimization by modifying the code above as follows: + + S2RegionTermIndexer::Options options; + options.set_index_contains_points_only(true); + S2RegionTermIndexer indexer(options); + +This reduces the number of query terms by approximately a factor of two, but +can only be used when an index contains only points. + +The complete example program looks like this: + +```c++ + // Create a set of "documents" to be indexed. Each document consists of a + // single point. + std::vector documents; + for (int docid = 0; docid < FLAGS_num_documents; ++docid) { + documents.push_back(S2Testing::RandomPoint()); + } + + // We use a hash map as our inverted index. The key is an index term, and + // the value is the set of "document ids" where this index term is present. + std::unordered_map> index; + + // Create an indexer suitable for an index that contains points only. + S2RegionTermIndexer::Options options; + options.set_index_contains_points_only(true); + S2RegionTermIndexer indexer(options); + static const char kPrefix[] = "s2:"; + + // Add the documents to the index. + for (int docid = 0; docid < documents.size(); ++docid) { + S2Point index_region = documents[docid]; + for (const auto& term : indexer.GetIndexTerms(index_region, kPrefix)) { + index[term].push_back(docid); + } + } + + // Convert the query radius to an angle representation. + S1Angle radius = + S1Angle::Radians(S2Earth::KmToRadians(FLAGS_query_radius_km)); + + // Count the number of documents (points) found in all queries. + int64 num_found = 0; + for (int i = 0; i < FLAGS_num_queries; ++i) { + // Choose a random center for query. + S2Cap query_region(S2Testing::RandomPoint(), radius); + + // Convert the query region to a set of terms, and compute the union of + // the document ids associated with those terms. + std::set candidates; + for (const auto& term : indexer.GetQueryTerms(query_region, kPrefix)) { + candidates.insert(index[term].begin(), index[term].end()); + } + + // "candidates" now contains all documents that intersect the query + // region, along with some documents that nearly intersect it. We can + // prune the results by retrieving the original "document" and checking + // the distance more precisely. + std::vector result; + for (int docid : candidates) { + if (!query_region.Contains(documents[docid])) continue; + result.push_back(docid); + } + // Now do something with the results (in this example we just count them). + num_found += result.size(); + } + printf("Found %lld points in %d queries\n", num_found, FLAGS_num_queries); +``` + +## Other S2Region Types + +As its name implies, `S2RegionTermIndexer` supports indexing and querying any +type of `S2Region`. In addition to points and discs (`S2Cap`), other useful +`S2Region` types include: + +* [`S2LatLngRect`](/devguide/basic_types#s2latlngrect): + a rectangle in latitude-longitude coordinates. +* [`S2Polyline`](/devguide/basic_types#s2polyline): a + polyline. +* [`S2Polygon`](/devguide/basic_types#s2polygon): a + polygon (can have holes and multiple shells). +* [`S2CellUnion`](/devguide/s2cell_hierarchy#s2cellunion): + a region approximated as a collection of `S2CellIds`. +* `S2ShapeIndexRegion` - an arbitrary collection of points, + polylines, and polygons. +* `S2ShapeIndexBufferedRegion` - like the above, but expanded by a + given radius. +* `S2RegionUnion` - the union of arbitrary other regions. +* `S2RegionIntersection` - the intersection of arbitrary other regions. + +For example, you could use `S2RegionTermIndexer` to index a set of polylines, +and then query which polylines intersect a given polygon. + +## Comparing Indexing Classes + +This tutorial gave a brief introduction to two useful indexing classes +(`S2PointIndex` and `S2RegionTermIndexer`). Here is a brief comparison of +their features: + +`S2PointIndex`: + +* Represents a collection of points in memory. +* Supports incremental updates. +* Can use `S2ClosestPointQuery` to find the points near a given target: + * Target types include points, polylines, polygons, and collections. + * Can also restrict results to a given `S2Region`. + +`S2RegionTermIndexer`: + +* A tool for adding spatial data to an information retrieval system. +* Can add one or more `S2Region` objects to each document. +* Queries return all documents that intersect a given `S2Region`. +* Requires an external system for storing and retrieving index terms. + +Another useful indexing class not described here is +[`S2ShapeIndex`](/devguide/s2shapeindex), which indexes an +arbitrary collection of points, polylines and polygons (collectively known as +*shapes*). `S2ShapeIndex` is the most generally useful of these classes for +working with geometry in memory. It compares to the classes above as follows: + +`S2ShapeIndex`: + +* Represents a collection of points, polylines, and polygons in memory. +* Supports incremental updates. +* Useful query classes include: + * `S2ContainsPointQuery`: returns the shape(s) that contain a given + point. + * [`S2ClosestEdgeQuery`](/devguide/s2closestedgequery): + returns the closest edges to a given point, polyline, polygon, or + geometry collection. + * `S2CrossingEdgeQuery`: returns the edge(s) that cross a given edge. + * `S2BooleanOperation`: computes boolean operations such as union, and + boolean predicates such as containment. + * `S2ShapeIndexRegion`: allows geometry collections to be approximated. + * `S2ShapeIndexBufferedRegion`: computes approximations that have been + expanded by a given radius. diff --git a/docs/devguide/degenerate_geometry.md b/docs/devguide/degenerate_geometry.md new file mode 100644 index 00000000..1f900656 --- /dev/null +++ b/docs/devguide/degenerate_geometry.md @@ -0,0 +1,584 @@ +--- +title: Degenerate Geometry in S2 +--- + +## Overview + +Unlike many libraries, S2 fully supports degenerate geometry. Degeneracies not +only do not cause errors, they have well-defined semantics including support for +open, semi-open, and closed boundaries. Supported degeneracy types include the +following: + + - *Point polylines* consisting of a single degenerate edge *AA* (where *A* + represents a vertex). + + - *Point loops* consisting of a single vertex *A*. Such loops may represent + either *point shells* or *point holes* according to whether the loop adds to + or subtracts from the surrounding region of the polygon. + + - *Sibling edge pairs* of the form \{*AB*, *BA*\}. Such sibling pairs may + represent either degenerate shells or degenerate holes according to whether + they add to or subtract from the surrounding region. The edges of a sibling + pair may belong to the same polygon loop (e.g. a loop *AB*) or to different + polygon loops (e.g. the polygon consisting of the loops *ABC* and *CBD*). + +Supporting such degeneracies has one big advantage, namely that geometry can be +simplified without completely losing small details. For example consider a +polygon representing a land region with many small lakes. Each lake would be +represented as a hole. When this polygon is simplified, some or even all of the +lakes may collapse to single points. Without support for degeneracies, such +lakes would disappear from the simplified result. If we allow degeneracies, on +the other hand, we can guarantee that for every point in the original region +there is a nearby corresponding point in the simplified region and vice versa. +In other words every lake in the input is still near some lake in the output. + +## Dimension Reduction: A Poor Substitute {#dimension-reduction} + +An alternative technique for handling degeneracies is *dimension reduction*, +whereby degenerate geometry of dimension 2 or 1 is replaced by non-degenerate +geometry of dimension 1 or 0. For example, the sibling pair \{*AB*, *BA*\} +would be replaced by a polyline edge (either *AB* or *BA*, chosen arbitrarily), +and similarly the degenerate polyline *AA* would be replaced by a single point +*A*. + +This technique has four very significant drawbacks compared to true degeneracy +support as implemented by the S2 library: + + - **Dimension reduction is only capable of handling** ***positive + degeneracies*** (e.g., polygon shells) as opposed to *negative degeneracies* + (e.g., polygon holes). Consider again the example of simplifying a land + region with many small lakes that collapse to single points. Since these + lakes are represented as holes in the polygon rather than shells, they cannot + be replaced by points but instead must be simply discarded. Similarly, + dimension reduction is not capable of representing negative degeneracies of + dimension 1 (e.g., a land polygon containing a river that when simplified has + its width collapse to zero). + + - **Dimension reduction changes the boundary of the affected geometry.** In + general the boundary of a polygon is defined to consist of its edges, the + boundary of a polyline is defined to consist of its two endpoints, and points + are defined to have no boundary at all. This implies that when the + degenerate sibling pair \{*AB*, *BA*\} is replaced by the polyline edge *AB*, + for example, the definition of its boundary changes from the two edges + \{*AB*, *BA*\} to the two points \{*A*, *B*\}. This is a very significant + difference, since it can change the result of virtually all the spatial + relationship predicates defined by OGC Simple Features model (which are based + on evaluating whether the interior, boundary, and exterior of one region + intersect the interior, boundary, and exterior of the other as summarized by + the `DE-9IM` matrix). + + - **Dimension reduction forces clients to deal with heterogeneous geometry as a + possible result of virtually any geometric operation.** For example, + intersecting two polygons may yield a collection of geometry consisting of + polygons, polylines, and points. Similarly, intersecting two polylines may + yield a collection of polylines and points. In contrast, the true degeneracy + support offered by S2 means that any boolean operation on two polygons is + guaranteed to yield a polygon, and any boolean operation on two polylines is + guaranteed to yield a polyline. This is not only simpler for clients but + also reduces the opportunities for bugs (since the cases resulting in + mixed-dimension geometry may be unexpected and rare). + + - **Dimension reduction makes it hard for clients to say how degeneracies + should be handled.** Some clients may wish to discard degeneracies, generate + errors, or convert degeneracies back into non-degenerate forms (e.g. + replacing point shells with tiny triangles). True degeneracy support makes + it easy to implement any of these models, and even to delegate this decision + to the functions that convert geometry to non-degeneracy-supporting formats + (such as geoJSON). Dimension reduction, on the other hand, makes it + difficult to distinguish degeneracies from legitimate lower-dimensional + objects. For example, imagine intersecting two collections of polylines and + polygons. Some output polylines will correspond to portions of input + polylines, while others may correspond to intersecting polygons that abut + along an edge. + +Other alternatives, such as extracting small details into separate objects to +protect them from simplification and then merging them back in afterwards, can +best be described as hacks. Not only are such algorithms virtually impossible +to make robust, they cannot offer the mathematical guarantees that can be +achieved by simplifying with a uniform error tolerance using degeneracies (see +below). + +## Degeneracies and Hausdorff Distance + +The fact that S2 supports both positive and negative degeneracies allows us to +make a very powerful guarantee when simplifying geometry. Recall that the +*Hausdorff distance* between two geometries *A* and *B* is defined as + +$$H(A, B) = {\rm max} \{ h(A, B), h(B, A) \}$$ + +where $$h(A, B)$$ is the *directed Hausdorff distance* + +$$h(A, B) = {\rm max}_{a \in A} \big\{ {\rm min}_{b \in B} \{ d(a, b) \} \big\} +. $$ + +**The Hausdorff distance is a measure of how different two geometries are;** in +particular it is an upper bound on how far away a point in one geometry can be +from the other geometry. + +The most important property of Hausdorff distance is that it can be used to +bound the distance measurement errors due to simplifying geometry. Suppose that +geometry *B* is a simplified version of geometry *A*, and that we wish to +measure the distance to some third geometry *X*. Then the error due to +simplification, i.e. the difference between $$d(A, X)$$ and $$d(B, X)$$, is +bounded by the Hausdorff distance between *A* and *B*: + +$$| d(A, X) - d(B, X) | \leq H(A, B) \, .$$ + +We can now state the guarantee that the S2 library makes when simplifying +geometry. Let geometry *A* be an arbitrary collection of points, polylines, and +polygons, and suppose that it is simplified using `S2Builder` with a snap radius +of *r* to yield a new geometry collection *B*. Then if degeneracies are allowed +when constructing *B*, we can guarantee that + +$$H(A, B) \leq r_{\rm e} \, .$$ + +The quantity $$r_{\rm e}$$ is called the *maximum edge deviation*, which +`S2Builder` ensures is at most 10% greater than the given snap radius *r* (i.e., +$$r_{\rm e} \leq 1.1 r$$). **In other words, the maximum possible distance +measurement error due to simplification is only slightly greater than the snap +radius used during simplification.** (Note that the additional 10% is a +guaranteed bound rather than a heuristic; it is only necessary because S2 works +with spherical geometry and would not be needed in a planar version of the +algorithm.) + +In fact, under the view that positive and negative degeneracies have +infinitesimal areas, the bound above applies not only to the whole geometries *A* +and *B*, but also to the interiors, boundaries, and exteriors of these +geometries: + +$$H({\rm Int}\, A, {\rm Int}\, B) \leq r_{\rm e}$$ + +$$H(\partial A, \partial B) \leq r_{\rm e}$$ + +$$H({\rm Ext}\, A, {\rm Ext}\, B) \leq r_{\rm e} \, .$$ + +These properties are the holy grail of simplification algorithms; it is a +mathematical way of saying that no details larger than the maximum edge devation +$$r_{\rm e}$$ have been lost. + +## Representation and Validity + +Degenerate geometry is not supported by the legacy classes `S2Polyline` and +`S2Polygon`. Instead the more modern, efficient, lightweight classes +`S2LaxPolylineShape` and `S2LaxPolygonShape` must be used. (The `Lax` in the +names refers to the fact that these classes allow degeneracies.) For example, +while `S2Polygon` requires that all loops have at least three vertices, +`S2LaxPolygonShape` supports loops with two, one, or even zero vertices. (The +zero-vertex case is called the *full loop* and represents the full sphere.) + +The `Lax` classes themselves do not have any validity requirements, but the +geometry they represent must satisfy certain conditions before it can be used in +S2 operations. Some operations are more restrictive than others; for example, +distance measurement has fewer requirements than `S2BooleanOperation`. You will +generally want to construct your geometry to satisfy the most restrictive rules, +outlined below. This can easily be accomplished by using `S2Builder` with an +appropriate output layer type (e.g. `S2LaxPolygonShape`). + +### Polylines + +Polylines must consist of either a single degenerate edge *AA*, or a sequence of +non-degenerate edges (e.g., *ABC*). Polylines with both degenerate and +non-degenerate edges are not allowed (e.g., *AABC* or *ABBC*). On the other +hand, polylines may self-intersect or have duplicate edges. + +Note that a point polyline is represented as two vertices (*AA*) whereas a point +polygon loop is represented as one vertex (*A*). Each of these objects thus +consists of a single degenerate edge AA, corresponding to the general rule that +an *n*-vertex loop has *n* edges, whereas an *n*-vertex polyline has *n-1* +edges. + +Here are a few examples of valid and invalid polylines: + + - *AB*: a valid polyline defining a single directed edge. + + - *ABCD*: a valid polyline consisting of three directed edges \{*AB*, *BC*, + *CD*\} + + - \[\] (the empty vertex sequence): a valid polyline defining no edges. + + - *AA*: a valid polyline defining a single degenerate edge at point *A*. + + - *ABBC*: invalid because degenerate and non-degenerate edges are used in the + same polyline. + + - *A*: invalid because it defines a vertex but no edges. (A degenerate edge at + *A* is represented as *AA*, while an empty polyline is represented as \[\].) + +### Polygons + +Polygons must consist of a set of oriented loops such that the interior of the +polygon is always on the left. Polygon edges may intersect only at their +endpoint vertices. A polygon may contain both an edge and its reverse edge +(i.e., *AB* and *BA*), but duplicate edges are not allowed (*AB* and *AB*). +Loops also must not contain repeated adjacent vertices (e.g., *ABBC*). + +Point loops consisting of a single vertex (*A*) are allowed, but such loops +cannot be incident to any other edges. The *full loop* (i.e., the loop with +zero vertices mentioned above) is allowed only if all the remaining loops (taken +together) define only degeneracies. + +Here are a few examples of valid and invalid polygons: + + - {} (the empty set of loops): a valid polygon containing no points (the empty + polygon). + + - \{*A*, *BC*\}: a valid polygon consisting of one point shell and one sibling + pair shell. + + - \{*AB*, *BAC*\}: invalid because the edge *BA* occurs twice (note that loop + *AB* consists of the two edges *AB* and *BA*). + + - \{*ABC*, *A*\}: invalid because the point shell *A* is incident to other + edges. + + - {full}: a valid polygon containing the entire sphere (the full polygon). + + - {full, *ABCB*\}: a valid polygon consisting of the entire sphere except for + two sibling pair holes (*AB* and *BC*). + + - {full, *ABC*\}: invalid because the full loop is used together with + non-degenerate geometry. + + - {full, *ABC*, *ACB*\}: a valid polygon consisting of the entire sphere except + for three sibling pair holes (*AB*, *BC*, and *CA*\}. Note that even though + the loops *ABC* and *ACB* are themselves non-degenerate, together they define + only degenerate geometry and therefore can be used with the full loop. + +### Mixed-Dimensional Geometry {#mixed-dimensional} + +The `S2ShapeIndex` class represents an arbitrary collection of points, +polylines, and polygons. Most S2 operations (e.g. `S2ClosestEdgeQuery`, +`S2RegionCoverer`, `S2ContainsPointQuery`) do not impose any additional validity +requirements on `S2ShapeIndex` beyond requiring that its component shapes are +valid. The notable exception to this rule is `S2BooleanOperation`, which +requires that **polygon interiors must be disjoint from all other geometry** +(including other polygon interiors). So for example, polygon interiors must not +intersect any point or polyline geometry. + +In order to support degeneracies, `S2BooleanOperation` requires in addition that +**duplicate polygon edges are not allowed.** Note that this rule is only +necessary when an `S2ShapeIndex` contains more than one polygon, since duplicate +edges are already disallowed within individual polygons. This additional rule +essentially requires that the collection of all polygons in each `S2ShapeIndex` +must satisfy the same validity requirements as a single polygon. (Also note +that more than one polygon is never needed, because polygons in S2 can have more +than one outer shell and therefore can represent an arbitrary subset of the +sphere.) + +Note that although duplicate polygon edges are not allowed, duplicate points and +polyline edges are permitted. So for example, an `S2ShapeIndex` containing the +points \{*A*, *A*\}, the polyline ABCABC, and the degenerate loop *AB* is a +valid input to `S2BooleanOperation`. This is true whether the boundary model is +open, semi-open, or closed. However an `S2ShapeIndex` input containing two +polygons each consisting of the degenerate loop *AB* would be invalid. + +Supporting points and polyline edges as multisets makes it much easier to +support time series such as GPS tracks. For example, a car driving around a +track might easily generate a polyline such as ABCABCABC. Similarly, this makes +it easier to reconstruct geometry that has been snapped and/or simplified, which +can cause distinct vertices to merge together. Note that if duplicate values +are not desired, they can easily be removed by choosing the appropriate +`S2Builder` output layer options. + +## Semantics of Degeneracies {#semantics} + +In general **degeneracies are modeled as infinitesimal regions of the +appropriate dimension.** So for example, a point shell *A* should be thought of +as a tiny closed loop in the vicinity of the vertex *A*, and a degenerate shell +loop *AB* (consisting of the two edges \{*AB*, *BA*\}) should be thought of as a +narrow sliver in the vicinity of edge *AB*. Similarly, the polyline *AA* should +be thought of as a tiny polyline in the vicinity of vertex *A*. We say "in the +vicinity of" because under certain boundary models, degeneracies may not contain +their nominal vertices and edges. For example, under the open boundary model +the point shell *A* does not contain the point *A*. Instead it should be +thought of as a tiny loop near the point *A*, so tiny that it does not contain +any representable points along its edges or in its interior. + +This model implies that **lower-dimensional degeneracies are infinitesimal +subsets of higher-dimensional degeneracies.** For example, the point *A* is +vanishingly small compared to the degenerate polyline *AA*. This is because the +point *A* is truly zero-dimensional, whereas the polyline *AA* represents a very +small one-dimensional set. As another example, consider a point *A*, a closed +polyline *AA*, and a closed point shell *A*. The intersection of the point *A* +and the closed polyline *AA* consists only of the point *A*, because it is an +infinitesimal subset of the polyline, and similarly the closed polyline *AA* +intersected with the closed polygon *A* yields only the polyline *AA*. + +It is also worth recalling that in the S2 library, **polyline and polygon edges +do not contain any interior points.** For example, the closed polygon shell +loop *AB* (consisting of the sibling edge pair \{*AB*, *BA*\}) does not contain +any points except *A* and *B*. This is because the S2 library models all +geometry as being infinitesimally perturbed such that edges do not pass through +any representable points on the sphere. This technique is known as *simulation +of simplicity* (Edelsbrunner and Muecke, 1990) and greatly reduces the number of +special cases that need to be handled when implementing robust geometric +algorithms. + +## Boundary Models + +The S2 library supports modeling polylines and polygons as having open, +semi-open, or closed boundaries. Here we briefly review these models and +discuss the treatment of degeneracies within them. + +Note that the polyline and polygon boundary models are specified independently. +The `PolylineModel` class specifies whether polylines contain their start and/or +end vertex, whereas the `PolygonModel` class specifies whether polygons contain +their vertices and edges. + +Also recall that the boundary model is specified for each operation rather than +being a property of the geometry itself. This allows more flexibility; for +example, if two geometries *A* and *B* intersect under the CLOSED boundary model +but not under the OPEN boundary model, this would mean that the boundaries of +the two geometries intersect but not their interiors. + +### Closed Boundary Model + +**A closed polyline contains all of its vertices.** For example, the closed +polyline *ABC* intersected with the three points \{*A*, *B*, *C*\} yields \{*A*, +*B*, *C*\}. + +Similarly, **a closed polygon contains all of its vertices, edges, and reversed +edges.** So for example, the closed polygon *ABC* intersected with the points +\{*A*, *B*, *C*\} yields \{*A*, *B*, *C*\}, and the closed polygon *ABC* +intersected with the six polylines \{*AB*, *BA*, *BC*, *CB*, *AC*, *CA*\} yields +\{*AB*, *BA*, *BC*, *CB*, *AC*, *CA*\}. + +These rules also apply to degeneracies. So for example, a closed point polyline +*AA* contains its vertex *A*, and a closed sibling pair shell *AB* contains its +vertices \{*A*, *B*\} and its edges \{*AB*, *BA*\}. + +Note that **certain degeneracies have no effect on point or edge containment; we +call such degeneracies *redundant*.** In the case of the closed boundary model, +**degenerate holes are redundant** because they do not exclude any points or +edges. For example, consider a polygon *ABC* containing a point hole *D*. +Since the boundary is closed, the polygon contains the point *D* whether or not +the hole is present. (The degenerate hole is still considered to exclude an +infinitesimal region from the polygon, it's simply that this region does not +contain any representable points.) + +Redundant degeneracies are considered to be part of the polygon's boundary and +therefore still affect distance measurement to the boundary. In fact, they may +be viewed as a method of expanding the polygon's boundary without affecting its +interior. They are treated in exactly the same way as non-degenerate regions in +boolean operations; for example, a redundant degeneracy intersected or unioned +with itself yields the original degeneracy, whereas a redundant degeneracy +subtracted from itself yields the empty set. + +### Open Boundary Model + +**An open polyline contain all of its vertices except the first and last.** For +example, the open polyline *ABC* intersected with the three points \{*A*, *B*, +*C*\} yields \{*B*\}. The only exception is if the +`polyline_loops_have_boundaries()` option is false and the first and last +vertices of a polyline are the same (e.g. *ABCDA*); in this case the first/last +vertex of the polyline loop is defined to be contained and its boundary is +defined to be empty. So for example, with this option the open polyline *ABCA* +intersected with \{*A*, *B*, *C*\} yields \{*A*, *B*, *C*\}. + +**An open polygon contains none of its vertices, edges, or reversed edges.** For +example, the open polygon *ABC* intersected with the points \{*A*, *B*, *C*\} +yields the empty set, and similarly the open polygon *ABC* intersected with the +six polylines \{*AB*, *BA*, *BC*, *CB*, *AC*, *CA*\} yields the empty set. + +These rules also apply to degeneracies. So for example, an open point polyline +*AA* contains no points at all, and an open sibling pair shell *AB* contains no +points or edges. In other words, just as degenerate polygon holes are redundant +in the closed boundary model, **degenerate polylines and degenerate polygon +shells are redundant in the open boundary model.** These types of degeneracies +do not affect point or edge containment, but simply add to the boundary of the +affected geometry. The geometry contains exactly the same set of points and +edges whether these degeneracies are present or not, however since the +degeneracies change the boundary they can still affect distance measurement. + +Also note that **geometric objects have exactly the same boundary under all +three boundary models.** So for example, the boundary of an open degenerate +polyline *AA* is the point *A*, even though the polyline is defined not to +contain that point. + +### Semi-Open Boundary Model + +The purpose of the semi-open boundary model is to allow geometry to **cover a +region without gaps or overlaps.** + +In order to achieve this, **a semi-open polyline contains all of its vertices +except the last.** For example, the semi-open polyline *ABC* intersected with +the points \{*A*, *B*, *C*\} yields \{*A*, *B*\}. The major advantage of this +model is that polylines can be split into pieces without affecting vertex +containment. For example, the single polyline *ABCDEFG* contains exactly the +same set of vertices as the three polylines \{*ABC*, *CDE*, *EFG*\}, and +furthermore each vertex is contained by exactly one of these polylines. + +Similarly, semi-open polygon point containment is defined such that **if several +semi-open polygons tile the region around a vertex, then exactly one of those +polygons contains that vertex.** This implies that a triangle *ABC* may or may +not contain each of its vertices \{*A*, *B*, *C*\}. However it ensures the +important property that **if a set of polygons tiles the entire S2 sphere, then +every point on the sphere is contained by exactly one polygon.** This property +can be very useful for certain algorithms. + +**Semi-open polygons contain all of their edges, but none of their reversed +edges.** For example, consider two abutting triangles *ABC* and *CBD*. The +edge *BC* is contained by *ABC* but not *CBD*, and the edge *CB* is contained +by *CBD* but not *ABC*. This property ensures **when a set of polygons tiles +the sphere, each polygon edge is contained by exactly one polygon.** + +This rule also ensures that if a polyline is intersected with a set of polygons +that tile the sphere, the resulting polylines can be unioned to obtain the +original polyline without gaps or overlaps. For example, let *ABCD* be a small +square polygon and consider its complement *DCBA*. These two polygons do not +intersect and their union is the entire S2 sphere. Now consider the zig-zag +polyline *ABDC*: its intersection with the polygon *ABCD* is the polyline *ABCD* +(, whereas its intersection with the polygon *DCBA* is *DC*. Together these two +polylines \{*ABD*, *DC*\} can be unioned to obtain the original polyline *ABDC* +without gaps or overlaps. + +**All types of degeneracies are redundant in the semi-open model.** In other +words, degeneracies never affect point or edge containment but simply add to the +boundary of the geometry. For example, a semi-open point polyline *AA* does not +contain the point *A* (just like the open model), semi-open degenerate polygon +shells do not contain any points or edges (just like the open model), and +semi-open degenerate polygon holes do not exclude any points or edges (just like +the closed model). All of these degeneracies simply add to the boundary of the +geometry without affecting its interior or exterior. The semi-open model is +purest of the three models in that shells and holes are treated symmetrically. + +## Boolean Operations + +Many properties of boolean operations involving degeneracies follow naturally +from the definitions above. Here we summarize the additional principles that +define the results of such operations (intersection, union, difference, and +symmetric difference). + +**Operations on geometry of a given dimension yields geometry of the same +dimension.** So for example, a polyline *AB* intersected with a polyline *BC* +under the closed model yields the degenerate polyline *BB* rather than a point. +Similarly a polygon *ABC* intersected with an abutting polygon *CBD* under the +closed model yields the degenerate polygon shell *BC*, not a polyline. And +symmetrically, a polygon hole *ACB* unioned with an abutting polygon hole *BCD* +under the open model yields a degenerate polygon hole *BC* (which has no +alternative representation). This behavior is exactly what allows S2 to avoid +the serious problems of [dimension reduction](#dimension-reduction). + +**In all boundary models, subtracting an object from itself yields the empty +set.** This is true even for redundant degeneracies that contain no points or +edges, such as subtracting a point polyline *AA* from itself in the open +boundary model. + +**In all boundary models, intersecting or unioning an object with itself yields +the original object(s).** The reason for the plural is that points and polyline +edges are treated as multisets. For example, if a point *A* is intersected with +another point *A*, the result is two points \{*A*, *A*\}. (Note that such +duplicates can easily be filtered away by choosing appropriate options in the +`S2Builder` output layer.) Again, this is true even for redundant degeneracies +that contain no points or edges, so for example a point shell *A* intersected +with itself in the open boundary model yields the same point shell *A*, and a +point polyline *AA* intesected with itself in the open boundary model yields two +point polylines \{*AA*, *AA*\} (neither of which contains the point *A*). + +**Intersecting geometry of different dimensions yields only the +lower-dimensional geometry.** This is consistent with the principles described +under [Semantics of Degeneracies](#semantics) above. For example, intersecting +a point *A* with a closed polygon *ABC* yields only the point *A*, not a point +polygon shell *A* as well. Of course, if the input geometries do not intersect +then the result is empty. For example intersecting a point *A* with an open +polygon *ABC* yields nothing at all. + +**Unioning geometry of different dimensions yields only the higher-dimensional +geometry.** So for example, unioning a point *A* with a closed polyline *AA* +yields only the polyline *AA*. Of course if the input geometries do not +intersect then the result is simply the collection of both inputs. For example +the union of a point *A* with an open or semi-open polyline *AA* is the +collection of both objects \{*A*, *AA*\} since they do not intersect. + +**Unioning points with other points, or polylines with other polylines, yields a +multiset.** For example, the union of two identical points *A* is the multiset +\{*A*, *A*\}, and the union of two identical polylines *ABC* is the multiset +\{*ABC*, *ABC*\}. As [noted earlier](#mixed-dimensional) this behavior makes it +possible to reconstruct time series such as GPS tracks accurately. Unwanted +duplicates can be filtered by choosing appropriate options in the `S2Builder` +output layer. + +### Semi-Open Semantics + +The semi-open boundary model has a few properties that deserve further +explanation. + +**In the semi-open model, complementary degeneracies do not intersect.** So for +example, a point shell *A* does not intersect a point hole *A*, and a sibling +pair shell *AB* does not intersect a sibling pair hole *AB*. + +**In the semi-open model, a point polyline *AA* or point shell *A* intersects +another different geometry if and only if that other geometry contains the point +*A*.** So for example, a point shell *A* intersects a semi-open polygon *ABC* if +and only if *ABC* contains its vertex *A* under the usual semi-open rules. In +other words, positive point degeneracies are treated exactly like points (except +when a degeneracy is intersected with itself---see above). + +The purpose of this rule is to ensure that the defining property of the +semi-open model (i.e. the ability to cover a region without overlaps or gaps) +applies to point degeneracies as well as points. For example, if a degenerate +polyline *AA* or degenerate shell *A* is intersected with a set of semi-open +polygons that tile the region around a shared vertex *A*, then exactly one of +those polygons contains that degeneracy. Similarly, if a semi-open polyline +*ABC* is expressed as the disjoint union of sub-polylines \{*AB*, *BC*\} then +exactly one of those polylines contains the degenerate polyline *BB* (in +particular, *BC* contains *BB* while *AB* does not). + +As noted above, the only way that point degeneracies are handled differently + than points is with respect to self-intersection. Specifically, in the + semi-open model: + + - A point *A* intersected with a polyline *AA* or point shell *A* is empty + (since the latter objects do not contain any points). + - A polyline *AA* intersected with a point shell *A* is empty (since the point + shell contains no points or polylines). + - And yet, a point *A*, polyline *AA*, or point shell *A* intersected with + itself yields the original object. + +These statements are not contradictory; there is in fact a valid representation +of these objects as infinitesimal regions that yields this behavior. + +The complementary property is that **the union of a point hole *A* with another +different geometry eliminates the hole if and only if that other geometry +contains the point *A*.** So for example, a point hole *A* unioned with the +semi-open polygon *ABC* eliminates the hole from the result if and only if *ABC* +contains its vertex *A*. + +A similar property for edge degeneracies is that **sibling pair shells do not +intersect abutting polygons.** So for example, a shell *BC* does not intersect +the polygon *ABC* or *CBD*. The best way to think about this is that the +intersection region is only one-half of the sibling pair degeneracy region and +therefore cannot be represented. The complementary property is that **the union +of a sibling pair hole with an abutting polygon eliminates the hole**. So for +example, the union of the polygon *ABC* with the full sphere except for the hole +*BC* is the full sphere. + +### Limitations + +It is important to realize that in some cases the true mathematical result of an +operation cannot be represented. This is true even when no degeneracies are +present. + +For example, consider a triangle *ABC* that contains a nested triangle *DEF*, +and suppose that *DEF* is subtracted from *ABC* in the CLOSED boundary model. +Ideally the result would contain the edges \{*AB*, *BC*, *CA*\} but not the +edges \{*DE*, *EF*, *FD*\} (since the latter three edges belong to *DEF*, which +has been subtracted). However it is not possible to represent this in the +CLOSED model since all of these edges are on the boundary of the result. + +Similarly, results involving degeneracies sometimes cannot be exactly +represented. In such cases the preference is always to correctly represent the +**interior** of the result (including infinitesimal interiors) even when the +boundary of the result cannot be perfectly represented. + +For example, consider subtracting a point shell *D* from the interior of a +closed polygon *ABC*. The result is defined to be the polygon \{*ABC*, *D*\} +where *D* is now a point hole. This may seem strange since (1) the result +still contains the point *D* because the boundary is closed, and (2) there is +no obvious reason to include *D* in the result because point holes exclude no +points under the closed model (i.e., they are redundant). And yet this is the +**only** correct result under the S2 degeneracy model, because the original +point shell *D* is considered to include an infinitesimal two-dimensional area +near the point *D* in addition to *D* itself. The only way to exclude this +infinitesimal two-dimensional area from the result is to include *D* as a point +hole. diff --git a/docs/devguide/examples/coverings.md b/docs/devguide/examples/coverings.md new file mode 100644 index 00000000..73e2ace2 --- /dev/null +++ b/docs/devguide/examples/coverings.md @@ -0,0 +1,152 @@ +--- +title: S2 Covering Examples +--- + +Here are some examples of the "cell coverage" feature of the new +spherical geometry package. We pick some point on the sphere, and create +a circular cap of some size around that point. We then cover the cap +with some number of "cells". The result is then projected onto a plane +tangent to the center of the cap and parallel to the surface of the +sphere. The cells are drawn in blue and the cap itself is drawn in red. + +This point happens to be centered on the Google offices in Kirkland, +with a 5 km radius cap. When you try to cover that cap with just one +cell, you often end up with a large amount of extra space. Notice the +large cell compared to the circle. + +![](img/kirkland_1.gif) + +Same configuration, except we allow the use of 2 cells. + +![](img/kirkland_2.gif) + +Same configuration, except we allow the use of 3 cells. One of them is +quite small. + +![](img/kirkland_3.gif) + +Same configuration, except we allow the use of 4 cells. I've changed the +scale factor so we can see more detail. + +![](img/kirkland_4.gif) + +Same configuration, except we allow the use of 5 cells. + +![](img/kirkland_5.gif) + +Same configuration, except we allow the use of 6 cells. + +![](img/kirkland_6.gif) + +Same configuration, except we allow the use of 10 cells. + +![](img/kirkland_10.gif) + +Same configuration, except we allow the use of 20 cells. + +![](img/kirkland_20.gif) + +Same configuration, except we allow the use of 50 cells. + +![](img/kirkland_50.gif) + +Same configuration, except we allow the use of 500 cells. The scale has +been changed again to show more detail. + +![](img/kirkland_500.gif) + +This is a point at an edge between two faces of the cubic tiling, with a +cap size of about 500 km and limit of 4 cells. The ratio of "area +covered by the cells" to "area covered by the cap" is about 10, which is +about the worst case value when you have at least 4 cells. (The ratio +can be quite large, on the order of 10\^15, if you limit the covering to +fewer than 4 cells.) + +![](img/edge_4.gif) + +Same configuration, except we allow the use of 6 cells. + +![](img/edge_6.gif) + +Same configuration, except we allow the use of 50 cells. + +![](img/edge_50.gif) + +Same configuration, except we allow the use of 500 cells. + +![](img/edge_500.gif) + +This is a point at the corner where three faces meet. The cap radius is +10 km, and we're using 10 cells. + +![](img/corner_10.gif) + +Same configuration, except we allow the use of 20 cells. + +![](img/corner_20.gif) + +Same configuration, except we allow the use of 50 cells. + +![](img/corner_50.gif) + +Same configuration, except we allow the use of 200 cells. + +![](img/corner_200.gif) + +Same configuration, except we allow the use of 1000 cells. Note that in +fact the code did not use all 1000 cells: this covering only uses 706 +cells. There was a point where the code had 1000 candidate cells, but +after pruning the code merged some adjacent sets of cells into a smaller +number of larger cells. The details of the covering did not change, it +is merely a more efficient representation. + +![](img/corner_1000.gif) + +This shows a latitude/longtitude rectangle (instead of a circular cap) +which is roughly placed over Washington state. (The coordinate system +that we're using for the project is not oriented to put North up, which +is why it looks a bit weird.) Here is a covering with 10 cells.. + +![](img/washington_10.gif) + +Same configuration, except we allow the use of 20 cells. + +![](img/washington_20.gif) + +Same configuration, except we allow the use of 100 cells. + +![](img/washington_100.gif) + +This shows a latitude/longtitude rectangle extends from 60 degrees north +to 80 degrees north, and from -170 degrees to +170 degrees longitude. +The covering is limited to 8 cells. Notice that the hole in the middle +is completely covered. + +![](img/polar_8.gif) + +Same configuration, except we allow the use of 20 cells. + +![](img/polar_20.gif) + +Same configuration, except we allow the use of 100 cells. Now we have +sufficient numbers of cells to model the hole, but not the gap near the +date line. + +![](img/polar_100.gif) + +Same configuration, except we allow the use of 500 cells. + +![](img/polar_500.gif) + +Finally, here are some examples of geographic coverings. This is a covering of +Hawaii using 25 cells. + +![](img/hawaii.gif) + +Here is a covering of Florida using 22 cells. + +![](img/florida1.gif) + +And Florida using 152 cells. + +![](img/florida2.gif) diff --git a/docs/devguide/examples/img/corner_10.gif b/docs/devguide/examples/img/corner_10.gif new file mode 100644 index 00000000..ea19c14e Binary files /dev/null and b/docs/devguide/examples/img/corner_10.gif differ diff --git a/docs/devguide/examples/img/corner_1000.gif b/docs/devguide/examples/img/corner_1000.gif new file mode 100644 index 00000000..23d4c64d Binary files /dev/null and b/docs/devguide/examples/img/corner_1000.gif differ diff --git a/docs/devguide/examples/img/corner_20.gif b/docs/devguide/examples/img/corner_20.gif new file mode 100644 index 00000000..f7127cd3 Binary files /dev/null and b/docs/devguide/examples/img/corner_20.gif differ diff --git a/docs/devguide/examples/img/corner_200.gif b/docs/devguide/examples/img/corner_200.gif new file mode 100644 index 00000000..f94a39da Binary files /dev/null and b/docs/devguide/examples/img/corner_200.gif differ diff --git a/docs/devguide/examples/img/corner_50.gif b/docs/devguide/examples/img/corner_50.gif new file mode 100644 index 00000000..ca6cb4ab Binary files /dev/null and b/docs/devguide/examples/img/corner_50.gif differ diff --git a/docs/devguide/examples/img/edge_4.gif b/docs/devguide/examples/img/edge_4.gif new file mode 100644 index 00000000..a6bbed50 Binary files /dev/null and b/docs/devguide/examples/img/edge_4.gif differ diff --git a/docs/devguide/examples/img/edge_50.gif b/docs/devguide/examples/img/edge_50.gif new file mode 100644 index 00000000..965a3a6f Binary files /dev/null and b/docs/devguide/examples/img/edge_50.gif differ diff --git a/docs/devguide/examples/img/edge_500.gif b/docs/devguide/examples/img/edge_500.gif new file mode 100644 index 00000000..1707798f Binary files /dev/null and b/docs/devguide/examples/img/edge_500.gif differ diff --git a/docs/devguide/examples/img/edge_6.gif b/docs/devguide/examples/img/edge_6.gif new file mode 100644 index 00000000..60696afa Binary files /dev/null and b/docs/devguide/examples/img/edge_6.gif differ diff --git a/docs/devguide/examples/img/florida1.gif b/docs/devguide/examples/img/florida1.gif new file mode 100644 index 00000000..9885e8ce Binary files /dev/null and b/docs/devguide/examples/img/florida1.gif differ diff --git a/docs/devguide/examples/img/florida2.gif b/docs/devguide/examples/img/florida2.gif new file mode 100644 index 00000000..d0820062 Binary files /dev/null and b/docs/devguide/examples/img/florida2.gif differ diff --git a/docs/devguide/examples/img/hawaii.gif b/docs/devguide/examples/img/hawaii.gif new file mode 100644 index 00000000..298182d5 Binary files /dev/null and b/docs/devguide/examples/img/hawaii.gif differ diff --git a/docs/devguide/examples/img/kirkland_1.gif b/docs/devguide/examples/img/kirkland_1.gif new file mode 100644 index 00000000..24e1811f Binary files /dev/null and b/docs/devguide/examples/img/kirkland_1.gif differ diff --git a/docs/devguide/examples/img/kirkland_10.gif b/docs/devguide/examples/img/kirkland_10.gif new file mode 100644 index 00000000..aaf10aa6 Binary files /dev/null and b/docs/devguide/examples/img/kirkland_10.gif differ diff --git a/docs/devguide/examples/img/kirkland_2.gif b/docs/devguide/examples/img/kirkland_2.gif new file mode 100644 index 00000000..629a31e5 Binary files /dev/null and b/docs/devguide/examples/img/kirkland_2.gif differ diff --git a/docs/devguide/examples/img/kirkland_20.gif b/docs/devguide/examples/img/kirkland_20.gif new file mode 100644 index 00000000..b0cb6663 Binary files /dev/null and b/docs/devguide/examples/img/kirkland_20.gif differ diff --git a/docs/devguide/examples/img/kirkland_3.gif b/docs/devguide/examples/img/kirkland_3.gif new file mode 100644 index 00000000..fa810824 Binary files /dev/null and b/docs/devguide/examples/img/kirkland_3.gif differ diff --git a/docs/devguide/examples/img/kirkland_4.gif b/docs/devguide/examples/img/kirkland_4.gif new file mode 100644 index 00000000..ff3d0dcd Binary files /dev/null and b/docs/devguide/examples/img/kirkland_4.gif differ diff --git a/docs/devguide/examples/img/kirkland_5.gif b/docs/devguide/examples/img/kirkland_5.gif new file mode 100644 index 00000000..e3e1ab3b Binary files /dev/null and b/docs/devguide/examples/img/kirkland_5.gif differ diff --git a/docs/devguide/examples/img/kirkland_50.gif b/docs/devguide/examples/img/kirkland_50.gif new file mode 100644 index 00000000..3ecb889b Binary files /dev/null and b/docs/devguide/examples/img/kirkland_50.gif differ diff --git a/docs/devguide/examples/img/kirkland_500.gif b/docs/devguide/examples/img/kirkland_500.gif new file mode 100644 index 00000000..11ca930a Binary files /dev/null and b/docs/devguide/examples/img/kirkland_500.gif differ diff --git a/docs/devguide/examples/img/kirkland_6.gif b/docs/devguide/examples/img/kirkland_6.gif new file mode 100644 index 00000000..e5a036c5 Binary files /dev/null and b/docs/devguide/examples/img/kirkland_6.gif differ diff --git a/docs/devguide/examples/img/polar_100.gif b/docs/devguide/examples/img/polar_100.gif new file mode 100644 index 00000000..eba23b05 Binary files /dev/null and b/docs/devguide/examples/img/polar_100.gif differ diff --git a/docs/devguide/examples/img/polar_20.gif b/docs/devguide/examples/img/polar_20.gif new file mode 100644 index 00000000..41411c1b Binary files /dev/null and b/docs/devguide/examples/img/polar_20.gif differ diff --git a/docs/devguide/examples/img/polar_500.gif b/docs/devguide/examples/img/polar_500.gif new file mode 100644 index 00000000..5431bccc Binary files /dev/null and b/docs/devguide/examples/img/polar_500.gif differ diff --git a/docs/devguide/examples/img/polar_8.gif b/docs/devguide/examples/img/polar_8.gif new file mode 100644 index 00000000..35d6881d Binary files /dev/null and b/docs/devguide/examples/img/polar_8.gif differ diff --git a/docs/devguide/examples/img/washington_10.gif b/docs/devguide/examples/img/washington_10.gif new file mode 100644 index 00000000..77794e9e Binary files /dev/null and b/docs/devguide/examples/img/washington_10.gif differ diff --git a/docs/devguide/examples/img/washington_100.gif b/docs/devguide/examples/img/washington_100.gif new file mode 100644 index 00000000..bcf0e19f Binary files /dev/null and b/docs/devguide/examples/img/washington_100.gif differ diff --git a/docs/devguide/examples/img/washington_20.gif b/docs/devguide/examples/img/washington_20.gif new file mode 100644 index 00000000..22aedaaa Binary files /dev/null and b/docs/devguide/examples/img/washington_20.gif differ diff --git a/docs/devguide/img/hilbert-cell-center.gif b/docs/devguide/img/hilbert-cell-center.gif new file mode 100644 index 00000000..f41007f3 Binary files /dev/null and b/docs/devguide/img/hilbert-cell-center.gif differ diff --git a/docs/devguide/img/hilbert-figure.gif b/docs/devguide/img/hilbert-figure.gif new file mode 100644 index 00000000..2fb83ad4 Binary files /dev/null and b/docs/devguide/img/hilbert-figure.gif differ diff --git a/docs/devguide/img/nesting-animation.gif b/docs/devguide/img/nesting-animation.gif new file mode 100644 index 00000000..09f87be1 Binary files /dev/null and b/docs/devguide/img/nesting-animation.gif differ diff --git a/docs/devguide/img/nesting-equator-rings0.png b/docs/devguide/img/nesting-equator-rings0.png new file mode 100644 index 00000000..bb5000af Binary files /dev/null and b/docs/devguide/img/nesting-equator-rings0.png differ diff --git a/docs/devguide/img/nesting-equator-rings1.png b/docs/devguide/img/nesting-equator-rings1.png new file mode 100644 index 00000000..0cbd003c Binary files /dev/null and b/docs/devguide/img/nesting-equator-rings1.png differ diff --git a/docs/devguide/img/s2cell_global.jpg b/docs/devguide/img/s2cell_global.jpg new file mode 100644 index 00000000..ff9e3e2f Binary files /dev/null and b/docs/devguide/img/s2cell_global.jpg differ diff --git a/docs/devguide/img/s2curve-large.gif b/docs/devguide/img/s2curve-large.gif new file mode 100644 index 00000000..26a67eba Binary files /dev/null and b/docs/devguide/img/s2curve-large.gif differ diff --git a/docs/devguide/img/s2curve-small.gif b/docs/devguide/img/s2curve-small.gif new file mode 100644 index 00000000..f75ca234 Binary files /dev/null and b/docs/devguide/img/s2curve-small.gif differ diff --git a/docs/devguide/img/s2curve-tiny.gif b/docs/devguide/img/s2curve-tiny.gif new file mode 100644 index 00000000..0c1fd7c8 Binary files /dev/null and b/docs/devguide/img/s2curve-tiny.gif differ diff --git a/docs/devguide/img/s2hierarchy.gif b/docs/devguide/img/s2hierarchy.gif new file mode 100644 index 00000000..711aa35a Binary files /dev/null and b/docs/devguide/img/s2hierarchy.gif differ diff --git a/docs/devguide/img/s2shapeindex_example.png b/docs/devguide/img/s2shapeindex_example.png new file mode 100644 index 00000000..b7cc8e31 Binary files /dev/null and b/docs/devguide/img/s2shapeindex_example.png differ diff --git a/docs/devguide/img/xyz_to_uv_to_st.gif b/docs/devguide/img/xyz_to_uv_to_st.gif new file mode 100644 index 00000000..7ec9fa45 Binary files /dev/null and b/docs/devguide/img/xyz_to_uv_to_st.gif differ diff --git a/docs/devguide/index.md b/docs/devguide/index.md new file mode 100644 index 00000000..ae0248b0 --- /dev/null +++ b/docs/devguide/index.md @@ -0,0 +1,26 @@ +--- +title: Developer Guide +--- + +Welcome to the S2 Developer Guide! This guide provides conceptual and +programming documentation for working with the S2 Geometry Library. + +If you haven't already done so, make sure you read our +[Conceptual Overview](/about/overview). + +S2 consists of source code repositories for C++ and Java. Upon +completing the overview, run through the Quickstart: + +* [C++ Quickstart](cpp/quickstart) + +(A Java Quickstart is envisioned soon.) + +Ready for More? You probably want to read our programming concepts: + +* Read about [S2 Cells](s2cell_hierarchy) +* Learn about the [Basic Types](basic_types) +* Understand the main [S2ShapeIndex](s2shapeindex) Class +* Discover [Finding Nearby Edges](s2closestedgequery) +* Determine how to [Interpolate Along Chains](s2chain_interpolation_query) +* Understand why [Degenerate Geometry](degenerate_geometry) is useful +* Read the [S2 C++ Style Guide](coding_style) diff --git a/docs/devguide/java/benchmarks.md b/docs/devguide/java/benchmarks.md new file mode 100644 index 00000000..16e134ce --- /dev/null +++ b/docs/devguide/java/benchmarks.md @@ -0,0 +1,293 @@ +--- +title: Java Benchmarks for S2 +--- + + + + + +## Overview + +The Java implementation of S2 has a collection of benchmarks that measure +performance of many important parts of the library. Their intent is to provide +accurate data to guide work on performance improvements, and also to detect +unwanted performance regressions. + + + +The benchmark code can be found in +google3/javatests/com/google/common/geometry/benchmarks . + + +## Benchmark Organization + +The Java S2 benchmarks are written using JMH, the Java Microbenchmark Harness +framework. Documentation on JMH is available from the [official site](https://github.com/openjdk/jmh). + +The S2 benchmarks are in multiple files, one for each S2 class. For example, +benchmarks for S2CellId are in +[S2CellIdBenchmark.java](https://source.corp.google.com/piper///depot/google3/javatests/com/google/common/geometry/benchmarks/S2CellIdBenchmark.java). +Each benchmark has a corresponding Blaze build target, with a +`Benchmark_jmh` suffix, e.g. `S2CellIdBenchmark_jmh`. + +Within each file, there are many individual JMH benchmarks, each of which has +the **@Benchmark** annotation. Most of these measures the performance of a +single method, while some may measure a sequence of operations on the class. +For example, the S2CellId Benchmark +[getAllNeighbors](https://source.corp.google.com/search?q=S2CellIdBenchmark%20getAllNeighbors&sq=package:piper%20file:%2F%2Fdepot%2Fgoogle3%20-file:google3%2Fexperimental) +measures the time to compute all neighbors of leaf cells. + +In a normal JMH run, the framework forks a new process which starts a new Java +virtual machine. Within that process, it calls each Benchmark method as many +times as possible until a target time has elapsed. Each call to the method is +called an **Invocation**. The whole sequence of calls until the target time +elapses is called a benchmark **Iteration**. + +JMH performs multiple Iterations for each Benchmark. The set of iterations is +called a **Trial**. When a trial completes, the JVM running it shuts down. +Multiple trials can be run for a benchmark (or all benchmarks) by increasing the +number of **Forks**. Every fork starts a new JVM and runs a complete Trial +within it. Statistics are computed across all Iterations in all Trials, or +forks, for the average time per invocation with a confidence interval. + +(If the number of forks is set to zero, everything runs in the same single JVM, +with one trial for each benchmark. This approach is not recommended for getting +accurate results, but it is faster, and is used in benchmark unit tests.) + +Before JMH runs the timed iterations in a Trial, it will run several **Warmup** +iterations. These run in the same way as the timed ones, but the results +are not used in computing the resulting statistics. Warmup iterations allow the +new JVM to fully optimize the code for maximum performance before the timed +iterations, which is important to get consistent, low-variance results. When you +look at the output of the run, you'll usually see something like the following +example, with the first few (warmup) iterations having somewhat lower +performance. This effect can be *very* dramatic for some code. + +``` +# Warmup Iteration 1: 34.030 ns/op +# Warmup Iteration 2: 29.849 ns/op +# Warmup Iteration 3: 27.266 ns/op +Iteration 1: 27.690 ns/op +Iteration 2: 27.213 ns/op +Iteration 3: 27.292 ns/op +Iteration 4: 27.351 ns/op +Iteration 5: 27.533 ns/op +``` + +If the first timed iteration of a benchmark is noticeably slower (a longer +average time per operation) than following ones, you should use more or longer +warmup iterations; otherwise, the final statistics will have lower accuracy and +a wider confidence interval. + +Within the code, Benchmark methods are usually contained within a **@State** +object. Benchmarks typically iterate through a variety of test data to provided +by the State to ensure the code is run with a representative range of values. + +The `State` also manages **@Parameters**. By default, JMH will run every +@Benchmark in a @State with every combination of @Parameters, and report the +results separately. Sometimes, you may not want to benchmark every possible +combination, as this can result in an explosion of cases. There are several +possible solutions to trim these number of cases. One such way is to specify +parameter values on the command line, as [described below](#commandline-parameters) + +However, if you *never* want to explore the entire parameter space, +you can define a Java enum for the combinations of parameters you want to +benchmark, where the enum instances each contain a set of values which are the +actual parameters. A simpler solution is to throw an exception in @Setup for +unwanted combinations of Parameters. This is effective if there is a constraint +like "Param A must always be less than Param B". + +## Running Benchmarks + +Running benchmarks on a local workstation or CloudTop is convenient and +reasonably accurate during development, but you should avoid running other +significant tasks such as Chrome or an IDE at the same time. Run benchmarks +for one class like: + +``` +export B=javatests/com/google/common/geometry/benchmarks +blaze run -c opt $B:S2LoopBenchmark_jmh +``` + +[Various flags](#flags) can be passed to JMH to control how the benchmark runs. + +For greater accuracy and repeatability, benchmarks can be run on dedicated +**Perflab** machines in Borg. To run on Perflab, you'll first need to complete +the quickstart instructions in the +[Perflab Tutorial](https://g3doc.corp.google.com/devtools/crosstool/perflab/g3doc/tutorial.md) +page. + +Then, you can run a benchmark on Perflab as follows: + +``` +export B=javatests/com/google/common/geometry/benchmarks +blaze run -c opt --run_under='perflab --constraints=arch=x86_64' \ + $B:S2CellIdBenchmark_jmh +``` + +When the job starts, URLs will be printed where you can track the progress of +the run on Borg, and find the final results in CNS when the job completes. +The output files have a TTL, which is currently 60 days. + +Note: you can ignore the message: `*Unable to locate a valid CitC path for CL +None, CitC workspace "None" and snapshot ID "None"*`. Under the covers, +`run_under` is using [mcd_perflab](go/perflab), which is designed for A/B +experiments. This message is just mcd_perflab saying it doesn't have a B side +for an A/B experiment.) + +Add the `--perfgate` flag if you want to upload benchmark results to PerfGate, +like this: + +``` +export B=javatests/com/google/common/geometry/benchmarks +blaze run -c opt --run_under='perflab --constraints=arch=x86_64' \ + $B:S2CellIdBenchmark_jmh -- --perfgate +``` + +## More Examples + +#### Run only one or more specific @Benchmarks for a class + +This technique is useful if you're working on improving or understanding +performance of a single S2 method or something similar, and just want to +repeatedly test the effect of changes you are making with the most applicable +benchmark. + +Add the name(s) of the @Benchmark methods to the end of the command line. The +strings you provide are actually matched as regular expressions to the fully +qualified method name. So, depending on the organization and naming of your +State and Benchmarks, you can select what you want in various ways. + +For example, these three commands all do the same thing: + +``` +export B=javatests/com/google/common/geometry/benchmarks + +blaze run -c opt $B:S2CellIdBenchmark_jmh com.google.common.geometry.benchmarks.S2CellIdBenchmark.EvenlySpacedLeafCellIdsState.toPointRaw +blaze run -c opt $B:S2CellIdBenchmark_jmh EvenlySpacedLeafCellIdsState.toPointRaw +blaze run -c opt $B:S2CellIdBenchmark_jmh toPointRaw +``` + +To run both the `toPointRaw` and `toPoint` benchmarks using a regular expression +match: + +``` +export B=javatests/com/google/common/geometry/benchmarks + +blaze run -c opt $B:S2CellIdBenchmark_jmh toPoint* +``` + +Providing just the name of a State object will match all the benchmarks within +it, like + +``` +export B=javatests/com/google/common/geometry/benchmarks +blaze run -c opt $B:S2CellIdBenchmark_jmh EvenlySpacedLeafCellIdsState +``` + +You can also specify multiple benchmark names: + +``` +export B=javatests/com/google/common/geometry/benchmarks +blaze run -c opt $B:S2EdgeUtilBenchmark_jmh robustCrossing edgeCrosser100RobustCrossings +``` + +#### Run a specific @Benchmark with specific parameter values{#commandline-parameters} + +This technique is useful (for example) if you are working to improve performance +for a specific scenario, perhaps large overlapping polygons, and don't need to +get results for other scenarios. + +Use the `-p` flag, repeatedly if needed, and with comma-separated values for +each parameter. This example sets values for two parameters: + +``` +blaze run -c opt $B:S2ClosestEdgeQueryBenchmark_jmh \ + FindClosestToSmallAbuttingIndex.findClosestEdge -- -p numIndexEdges=1024,65536 -p factory=FRACTAL_LOOP +``` + +## JMH Benchmark Flags{#flags} + +JMH has many command flags that can specify values for many of the settings +described above. It is convenient to provide defaults for these values using +annotations like @Warmup and @Measurement in the code, at the level of +individual @Benchmarks or @States, but those can be overridden from the command +line with the following flags: + +| Flag Example | Meaning | +| ---------------- | ----------------------------------------------------------| +| -wi=3 | Set the number of warmup iterations per trial. | +| -i=5 | Set the number of timed benchmark iterations per trial. | +| -f=2 | Set the number of forks, i.e. how many times JMH will fork +: : a new JVM process to run a set of benchmark trials. : +| -w=2s | Set the target time per warmup iteration to 2 seconds | +| -r=10s | Set the target time per benchmark iteration to 10 seconds | +| -to=10 | Set the timeout for iterations to 10 seconds. (note, no | +: : trailing 's') : +| -p arg=val1,val2 | Override the parameter values used for parameter "arg" to : +: : only the specific provided values. Can be repeated. : + +**Flags Example** + +Note that blaze differentiates flags for the blaze binary itself from flags for +the run target using a double-hyphen. + +To run the `S2ClosestEdgeQuery.FindClosestToSmallAbuttingIndex` benchmark with +five timed iterations, 10 seconds per measured iteration, the +`numIndexedEdges` @Param set to 64K, and other values using the values set in +the benchmark code: + +``` +blaze run -c opt $B:S2ClosestEdgeQueryBenchmark_jmh FindClosestToSmallAbuttingIndex -- -i=5 -r=10s -p numIndexEdges=65536 +``` + +## Developing Benchmarks + +**Test the benchmarks for a single class:** + +Testing a benchmark with "blaze test" invokes each @Benchmark one time, for +every combination of parameters. The blaze target names are like +`test_S2LoopBenchmark_jmh`. For example: + +``` +export B=javatests/com/google/common/geometry/benchmarks +blaze test -c opt $B:test_S2ClosestEdgeQueryBenchmark_jmh +``` + +"Blaze test" is supposed to be fast, but if some combinations of parameters take +a long time to run even once, test timeouts can be a problem. The "isUnitTest()" +method can be used in @State `setup()` methods to replace expensive parameter +values with fast ones. See the +[S2PolygonBenchmark](https://source.corp.google.com/piper///depot/google3/javatests/com/google/common/geometry/benchmarks/S2PolygonBenchmark.java) +for some examples. + +**Test all the benchmarks:** + +``` +export B=javatests/com/google/common/geometry/benchmarks +blaze test -c opt $B/... +``` + +## Comparing Java to C++ benchmarks + +Many (but not all) of the Java S2 benchmarks were designed to do exactly the +same operations as corresponding C++ benchmarks. This allows for comparing the +performance of the Java and C++ implementations. + +The C++ benchmarks are part of the unit tests. You can build one like: + +``` +blaze build -c opt --dynamic_mode=off --copt=-gmlt util/geometry:s2closest_edge_query_test +``` + +And run all the benchmarks in a unit test: + +``` +blaze-bin/util/geometry/s2closest_edge_query_test --benchmark_filter=all +``` diff --git a/docs/devguide/s2cell_hierarchy.md b/docs/devguide/s2cell_hierarchy.md new file mode 100644 index 00000000..ce9e23f4 --- /dev/null +++ b/docs/devguide/s2cell_hierarchy.md @@ -0,0 +1,1041 @@ +--- +title: S2 Cells +--- + +The S2 library defines a framework for decomposing the unit sphere into a +hierarchy of *cells*. Each cell is a quadrilateral bounded by four geodesics. +The top level of the hierarchy is obtained by projecting the six faces of a +cube onto the unit sphere, and lower levels are obtained by subdividing each +cell into four children recursively. For example, the following image shows +two of the six *face cells*, one of which has been subdivided several times: + +![](img/s2hierarchy.gif) + +Notice that the cell edges appear to be curved; this is because they are +*spherical geodesics*, i.e., straight lines on the sphere (similar to the +routes that airplanes fly). + +Each cell in the hierarchy has a *level*, defined as the number of times the +cell has been subdivided (starting with a face cell). Cells levels range from +0 to 30. The smallest cells at level 30 are called *leaf cells*; there are +6 * 430 of them in total, each about 1cm across on the Earth's +surface. (Details on the cell sizes at each level can be found on the [S2 +Cell Statistics page](/resources/s2cell_statistics).) + +The S2 hierarchy is useful for spatial indexing and for approximating regions +as a collection of cells. Cells can be used to represent both points and +regions: points are generally represented as leaf cells, while regions are +represented as collections of cells at any level(s). For example, here +is an approximation of Hawaii as a collection of 22 cells: + +![](examples/img/hawaii.gif) + +## S2CellId Numbering + +Each cell is uniquely identified by a 64-bit [`S2CellId`](#s2cellid). +The S2 cells are numbered in a special way in order to maximize +[locality of +reference](https://en.wikipedia.org/wiki/Locality_of_reference) when +they are used for spatial indexing (compared to other methods of +numbering the cells). + +In particular, the S2 cells are ordered sequentially along a +[space-filling curve](https://en.wikipedia.org/wiki/Space-filling_curve){:target="_blank"} +(a type of [fractal](https://en.wikipedia.org/wiki/Fractal){:target="_blank"}). +The particular curve used by S2 is called the *S2 space-filling curve*, +and consists of six Hilbert curves linked together to form a single +continuous loop over the entire sphere. Here is an illustration of the +S2 curve after 5 levels of subdivision: + +[![](img/s2curve-small.gif)](img/s2curve-large.gif) + +The yellow curve is an approximation of the S2 space-filling curve. It +is a single continuous loop with a fractal structure such that it +passes near every point on the sphere. (If you were to cut along the +yellow line, it would separate the sphere into two equal halves.) The +green lines show the boundaries of the S2 cells at levels 0, 1, 2, 3, +4, and 5, drawn as lines of different widths. (You can click on the +image to see a larger version.) + +The cells at level 5 are numbered in increasing order along this curve. +This leads to the property that *if the `S2CellIds` of two cells are +close together, then the cells are also close +together*.[^s2cell-locality] This property significantly improves locality +of reference when the cells are used for indexing. + +The remainder of this section gives further details about how the S2 cell +hierarchy is organized. (You don't need to understand this background +information in order to use the library -- you can [skip forward](#s2cellid) if +desired.) + +### Hilbert Curve + +The S2 curve is based on the Hilbert curve. The Hilbert curve is a function +from the unit interval [0,1] to the unit square [0,1]×[0,1] that is +*space-filling*, meaning that it visits every point in the unit square. The +Hilbert curve is continuous but not differentiable, and can be considered to +have infinite length. + +It is most easily defined as the limit of an iterative process that builds a +more detailed approximation of the curve at each step. Here is a part of the +first page of +[Hilbert's 1891 paper defining his curve](http://www.digizeitschriften.de/dms/img/?PPN=PPN235181684_0038&DMDID=dmdlog40){:target="_blank"} +(cf. +[Mark McClure](https://math.stackexchange.com/users/21361/mark-mcclure){:target="_blank"} +): + +![Figures from Hilbert's 1891 paper](img/hilbert-figure.gif) + +As you can see, the first iteration divides the unit square into 4 +smaller squares. The curve visits those squares in a particular order +that looks like an inverted "U" (figure 1). The second iteration takes +each square from the first iteration and divides it into 4 smaller +subsquares. The curve again visits those subsquares in a U-shaped +order, except that some of the U-shapes have been rotated and/or +reflected in order to link the curves together seamlessly (figure 2). +(The rotation/reflection rules are simple but we will not describe them +here.) Figure 3 shows the result after three iterations. + +How does this process define a mapping from the unit interval to the unit +square? Let + +> H : [0,1] → [0,1] × [0,1] + +be the Hilbert curve, and suppose that we want to evaluate H(s) for +the real number s = 0.5. We start by writing out the binary expansion of +s = 0.5: + +> s = 0.100000000... + +To find the corresponding point H(s) on the Hilbert curve, we group the +digits of "s" into pairs: + +> [10, 00, 00, 00, ...] + +Each pair of binary digits corresponds to a decimal number between 0 and 3, +and indicates which of the 4 subsquares to choose at each step of the +construction process. Note that unlike Hilbert's paper, we number the +subsquares of each square from 0 to 3 (in the order they are visited by the +curve) rather than 1 to 4. For example, in the case of s = 0.5, the first +group of digits "10" corresponds to the decimal number 2, which means that we +choose subsquare 2 during the first iteration (i.e., square 3 in Hilbert's +figure 1). The next group of digits is "00", meaning that we choose subsquare +0 during the second iteration (corresponding to square 9 in Hilbert's figure +2). Continuing in this way, for s = 0.5 the result looks something like this: + +![](img/hilbert-cell-center.gif) + +In the limit, this process converges to a point that is exactly at the +center of the square (the green dot in the figure). This implies that + +> H(0.5) = (0.5, 0.5) + +i.e., the midpoint of the Hilbert curve is at the exact center of unit +square.[^cell-center-parameters] The value of H(s) for any real number "s" can +be evaluated using a similar process. + +### S2 Space-Filling Curve + +The S2 space-filling curve is a function + +> S : [0,6] → S² + +i.e., from the interval [0,6] to the unit sphere. It is constructed by +mapping 6 copies of the Hilbert curve to the 6 faces of a unit cube, +reflecting and rotating the curves as necessary so that that they link +together seamlessly into a continuous loop. The cube is then mapped to +the unit sphere using a transformation that minimizes distortion. + +We can view the S2 curve parameter "s" as having the format: + +> s = [face].[face_pos] + +where "face" is a number in the range [0..5] that selects one of the +six cube faces, and "face_pos" is the Hilbert curve parameter on that +face. This mapping is defined such that the end of the curve on face i +exactly matches the start of the curve on face i+1 (including the +transition from face 5 to face 0). + +Here is an illustration (not fully accurate, but conceptually correct) +of the S2 curve on the unit cube (before projecting it to the sphere). +The cube has been unfolded and flattened, and shows the first 4 levels +of the S2 curve subdivision: + +![](img/s2cell_global.jpg) + +([Larger images of the individual faces can be found +here.](/resources/earthcube.md)) Note that the traversal +order of the odd-numbered faces is the mirror image of the even-numbered faces. + +When the cube is mapped onto the unit sphere, the result is the curve +shown earlier: + +![](img/s2curve-tiny.gif) + +### S2CellId Numbering (again) + +We now return to the question of how the S2 cells are numbered. Recall +that there are 31 levels of subdivision, ranging from 0 to 30. It +turns out that the cell numbering system at all levels follows a very +simple rule, namely: + +> The S2CellId of cell C is the S2 curve parameter at the center of C +> (scaled to obtain an integer). + +For example, consider the level 0 cell for cube face 2. The center of +this cell has the S2 curve parameter (in binary): + +> s = 010.100000000... + +where "010" is the binary representation of face 2 and ".1000..." is +the binary representation of 0.5 (i.e., the Hilbert curve parameter at +the center of this face, as discussed above). Scaling this by a factor +of 261, we obtain 0101000... in binary (where there are 60 +trailing zeros), or 0x5000000000000000 in hexadecimal. This is the +`S2CellId` for the level 0 cell on face 2. + +It turns out that all `S2CellIds` up to level 30 can be represented in this +way, i.e. the S2 curve parameters of the corresponding cell centers can be +represented exactly in 64 bits. Note that no two cells (even at different +levels) have the same center position, and therefore all the `S2CellIds` are +distinct. The cell center is sufficient to specify both the position **and** +the subdivision level of each cell. + +We can also look at this from a more practical point of view. The +`S2CellId` for a cell at level k always has the following structure: + +> s = [face] [child]k 1 060-2k + +where "face" is a 3-bit number (range 0..5) that selects a cube face, +and each of the k "child" values is a 2-bit number that selects one of +the 4 children during the subdivision process. The bit that follows +the last child is always 1, and all other bits are zero. For example: + + 010 10000...0 Face cell 2. + 001 10 100...0 Subcell 2 of face cell 1. + 100 11 01 1000...0 Subcell 1 of subcell 3 of face 4. + +This representation has the convenient property that the subdivision +level of a cell can easily be determined from the position of its +lowest-numbered 1 bit. + +### Coordinate Systems + +The process of mapping an `S2CellId` to a point on the unit sphere involves +several steps, each with its own coordinate system. + +* (cellid)
+ Cell id: A 64-bit encoding of a face and a Hilbert curve parameter on that + face, as discussed above. The Hilbert curve parameter implicitly encodes + both the position of a cell and its subdivision level. + +* (face, i, j)\ + Leaf-cell coordinates: The leaf cells are the subsquares that result after + 30 levels of Hilbert curve subdivision, consisting of a 230 + × 230 array on each face. "i" and "j" are integers in + the range [0, 230-1] that identify a particular leaf cell. The + (i, j) coordinate system is right-handed on every face, and the faces are + oriented such that Hilbert curves connect continuously from one face to + the next. + +* (face, s, t)\ + Cell-space coordinates: "s" and "t" are real numbers in the range [0,1] + that identify a point on the given face. For example, the point (s, t) = + (0.5, 0.5) corresponds to the center of the cell at level 0. Cells in + (s, t)-coordinates are perfectly square and subdivided around their center + point, just like the Hilbert curve construction. + +* (face, u, v)\ + Cube-space coordinates: To make the cells at each level more uniform in + size after they are projected onto the sphere, we apply a nonlinear + transformation of the form u=f(s), v=f(t) before projecting points onto + the sphere. This function also scales the (u,v)-coordinates so that each + face covers the *biunit square* [-1,1]×[-1,1]. Cells in + (u,v)-coordinates are rectangular, and are not necessarily subdivided + around their center point (because of the nonlinear transformation "f"). + +* (x, y, z)\ + Spherical point: The final `S2Point` is obtained by projecting the (face, + u, v) coordinates onto the unit sphere. Cells in (x,y,z)-coordinates are + quadrilaterals bounded by four spherical geodesic edges. + +The purpose of the nonlinear (face, s, t) → (face, u, v) transformation +is to make the cells roughly the same size on the sphere. This can be +visualized as follows: + +[![](img/xyz_to_uv_to_st.gif){height="400"}](img/xyz_to_uv_to_st.gif) + +This diagram shows a one-dimensional slice of cube face. Starting from the +right-hand size, an s- or t-coordinate in the range [0,1] is transformed +non-linearly into a u- or v-coordinate in the range [-1,1]. This point is +then projected onto the unit sphere. + +The purpose of the (face, s, t)→(face, u, v) transformation is +illustrated by the two angles marked in red. Although both angles are the +same size, and therefore correspond to the same distance on the unit sphere, +they have very different sizes when projected back to (u,v)-coordinates (i.e., +the top interval on the middle (u,v)-slice is much larger). The non-linear +transformation between (s,t) and (u,v) coordinates helps to correct this +situation (i.e., the intervals on the right-hand (s,t)-slice are closer in +size). + +Note that the (i, j), (s, t), and (u, v) coordinate systems are right-handed +on all six faces. + +## S2CellId + +The `S2CellId` class is a thin wrapper over +[the 64-bit S2CellId number](#s2cellid-numbering) that provides methods for +navigating the cell hierarchy (finding parents, children, containment tests, +etc). Since leaf cells are often used to represent points on the unit sphere, +the `S2CellId` class also provides methods for converting directly to and from +an `S2Point`. Here are its methods: + +```c++ +class S2CellId { + public: + static const int kFaceBits = 3; + static const int kNumFaces = 6; + static const int kMaxLevel = 30; // Valid levels: 0..kMaxLevel + static const int kMaxSize = 1 << kMaxLevel; + + // Although only 60 bits are needed to represent the index of a leaf + // cell, we need an extra bit in order to represent the position of + // the center of the leaf cell along the Hilbert curve. + static const int kPosBits = 2 * kMaxLevel + 1; + + explicit S2CellId(uint64 id); + + // The default constructor returns an invalid cell id. + S2CellId(); + static S2CellId None(); + + // Return an invalid cell id guaranteed to be larger than any + // valid cell id. Useful for creating indexes. + static S2CellId Sentinel(); + + // Return the cell corresponding to a given S2 cube face. + static S2CellId FromFace(int face); + + // Return a cell given its face (range 0..5), Hilbert curve position within + // that face (an unsigned integer with S2CellId::kPosBits bits), and level + // (range 0..kMaxLevel). The given position will be modified to correspond + // to the Hilbert curve position at the center of the returned cell. This + // is a static function rather than a constructor in order to indicate what + // the arguments represent. + static S2CellId FromFacePosLevel(int face, uint64 pos, int level); + + // Construct a leaf cell containing the given point "p". Usually there is + // exactly one such cell, but for points along the edge of a cell, any + // adjacent cell may be (deterministically) chosen. This is because + // S2CellIds are considered to be closed sets. The returned cell will + // always contain the given point, i.e. + // + // S2Cell(S2CellId(p)).Contains(p) + // + // is always true. The point "p" does not need to be normalized. + explicit S2CellId(const S2Point& p); + + // Construct a leaf cell containing the given normalized S2LatLng. + explicit S2CellId(const S2LatLng& ll); + + // Return the direction vector corresponding to the center of the given + // cell. The vector returned by ToPointRaw is not necessarily unit length. + // This method returns the same result as S2Cell::GetCenter(). + S2Point ToPoint() const; + S2Point ToPointRaw() const; + + // Return the S2LatLng corresponding to the center of the given cell. + S2LatLng ToLatLng() const; + + // The 64-bit unique identifier for this cell. + uint64 id() const; + + // Return true if id() represents a valid cell. + bool is_valid() const; + + // Which cube face this cell belongs to, in the range 0..5. + int face() const; + + // The position of the cell center along the Hilbert curve over this face, + // in the range 0..(2^kPosBits-1). + uint64 pos() const; + + // Return the subdivision level of the cell (range 0..kMaxLevel). + int level() const; + + // Return true if this is a leaf cell (more efficient than checking + // whether level() == kMaxLevel). + bool is_leaf() const; + + // Return true if this is a top-level face cell (more efficient than + // checking whether level() == 0). + bool is_face() const; + + // Return the child position (0..3) of this cell within its parent. + // REQUIRES: level() >= 1. + int child_position() const; + + // Return the child position (0..3) of this cell's ancestor at the given + // level within its parent. For example, child_position(1) returns the + // position of this cell's level-1 ancestor within its top-level face cell. + // REQUIRES: 1 <= level <= this->level(). + int child_position(int level) const; + + // Methods that return the range of cell ids that are contained + // within this cell (including itself). The range is *inclusive* + // (i.e. test using >= and <=) and the return values of both + // methods are valid leaf cell ids. + // + // These methods should not be used for iteration. If you want to + // iterate through all the leaf cells, call child_begin(kMaxLevel) and + // child_end(kMaxLevel) instead. Also see maximum_tile(), which can be used + // to iterate through cell ranges using cells at different levels. + // + // It would in fact be error-prone to define a range_end() method, because + // this method would need to return (range_max().id() + 1) which is not + // always a valid cell id. This also means that iterators would need to be + // tested using "<" rather that the usual "!=". + S2CellId range_min() const; + S2CellId range_max() const; + + // Return true if the given cell is contained within this one. + bool contains(S2CellId other) const; + + // Return true if the given cell intersects this one. + bool intersects(S2CellId other) const; + + // Return the cell at the previous level or at the given level (which must + // be less than or equal to the current level). + S2CellId parent() const; + S2CellId parent(int level) const; + + // Return the immediate child of this cell at the given traversal order + // position (in the range 0 to 3). This cell must not be a leaf cell. + S2CellId child(int position) const; + + // Iterator-style methods for traversing the immediate children of a cell or + // all of the children at a given level (greater than or equal to the current + // level). Note that the end value is exclusive, just like standard STL + // iterators, and may not even be a valid cell id. You should iterate using + // code like this: + // + // for(S2CellId c = id.child_begin(); c != id.child_end(); c = c.next()) + // ... + // + // The convention for advancing the iterator is "c = c.next()" rather + // than "++c" to avoid possible confusion with incrementing the + // underlying 64-bit cell id. + S2CellId child_begin() const; + S2CellId child_begin(int level) const; + S2CellId child_end() const; + S2CellId child_end(int level) const; + + // Return the next/previous cell at the same level along the Hilbert curve. + // Works correctly when advancing from one face to the next, but + // does *not* wrap around from the last face to the first or vice versa. + S2CellId next() const; + S2CellId prev() const; + + // This method advances or retreats the indicated number of steps along the + // Hilbert curve at the current level, and returns the new position. The + // position is never advanced past End() or before Begin(). + S2CellId advance(int64 steps) const; + + // Return the level of the "lowest common ancestor" of this cell and + // "other". Note that because of the way that cell levels are numbered, + // this is actually the *highest* level of any shared ancestor. Return -1 + // if the two cells do not have any common ancestor (i.e., they are from + // different faces). + int GetCommonAncestorLevel(S2CellId other) const; + + // Iterator-style methods for traversing all the cells along the Hilbert + // curve at a given level (across all 6 faces of the cube). Note that the + // end value is exclusive (just like standard STL iterators), and is not a + // valid cell id. + static S2CellId Begin(int level); + static S2CellId End(int level); + + // Methods to encode and decode cell ids to compact text strings suitable + // for display or indexing. Cells at lower levels (i.e. larger cells) are + // encoded into fewer characters. The maximum token length is 16. + // + // ToToken() returns a string by value for convenience; the compiler + // does this without intermediate copying in most cases. + // + // These methods guarantee that FromToken(ToToken(x)) == x even when + // "x" is an invalid cell id. All tokens are alphanumeric strings. + // FromToken() returns S2CellId::None() for malformed inputs. + string ToToken() const; + static S2CellId FromToken(const char* token, size_t length); + static S2CellId FromToken(const string& token); + + // Creates a debug human readable string. Used for << and available for direct + // usage as well. + string ToString() const; + + // Return the four cells that are adjacent across the cell's four edges. + // Neighbors are returned in the order defined by S2Cell::GetEdge. All + // neighbors are guaranteed to be distinct. + void GetEdgeNeighbors(S2CellId neighbors[4]) const; + + // Return the neighbors of closest vertex to this cell at the given level, + // by appending them to "output". Normally there are four neighbors, but + // the closest vertex may only have three neighbors if it is one of the 8 + // cube vertices. + // + // Requires: level < this->level(), so that we can determine which vertex is + // closest (in particular, level == kMaxLevel is not allowed). + void AppendVertexNeighbors(int level, std::vector* output) const; + + // Append all neighbors of this cell at the given level to "output". Two + // cells X and Y are neighbors if their boundaries intersect but their + // interiors do not. In particular, two cells that intersect at a single + // point are neighbors. + // + // Requires: nbr_level >= this->level(). Note that for cells adjacent to a + // face vertex, the same neighbor may be appended more than once. + void AppendAllNeighbors(int nbr_level, std::vector* output) const; +}; +``` + +See `s2cellid.h` for additional methods. + +## S2Cell + +An `S2Cell` is an `S2Region` object that represents a cell. Unlike `S2CellId`, +it views a cell as a representing a spherical quadrilateral rather than a point, +and it supports efficient containment and intersection tests. However, it is +also a more expensive representation (currently 48 bytes rather than 8). + +Here are its methods: + +```c++ +class S2Cell final : public S2Region { + public: + // The default constructor is required in order to use freelists. + // Cells should otherwise always be constructed explicitly. + S2Cell(); + + // An S2Cell always corresponds to a particular S2CellId. The other + // constructors are just convenience methods. + explicit S2Cell(S2CellId id); + + // Convenience constructors. The S2LatLng must be normalized. + explicit S2Cell(const S2Point& p); + explicit S2Cell(const S2LatLng& ll); + + // Return the cell corresponding to the given S2 cube face. + static S2Cell FromFace(int face); + + // Return a cell given its face (range 0..5), Hilbert curve position within + // that face (an unsigned integer with S2CellId::kPosBits bits), and level + // (range 0..kMaxLevel). The given position will be modified to correspond + // to the Hilbert curve position at the center of the returned cell. This + // is a static function rather than a constructor in order to indicate what + // the arguments represent. + static S2Cell FromFacePosLevel(int face, uint64 pos, int level); + + S2CellId id() const; + int face() const; + int level() const; + int orientation() const; + bool is_leaf() const; + + // Return the k-th vertex of the cell (k = 0,1,2,3). Vertices are returned + // in CCW order (lower left, lower right, upper right, upper left in the UV + // plane). The points returned by GetVertexRaw are not normalized. + S2Point GetVertex(int k) const; + S2Point GetVertexRaw(int k) const; + + // Return the inward-facing normal of the great circle passing through + // the edge from vertex k to vertex k+1 (mod 4). The normals returned + // by GetEdgeRaw are not necessarily unit length. + S2Point GetEdge(int k) const; + S2Point GetEdgeRaw(int k) const; + + // If this is not a leaf cell, set children[0..3] to the four children of + // this cell (in traversal order) and return true. Otherwise returns false. + // This method is equivalent to the following: + // + // for (pos=0, id=child_begin(); id != child_end(); id = id.next(), ++pos) + // children[pos] = S2Cell(id); + // + // except that it is more than two times faster. + bool Subdivide(S2Cell children[4]) const; + + // Return the direction vector corresponding to the center in (s,t)-space of + // the given cell. This is the point at which the cell is divided into four + // subcells; it is not necessarily the centroid of the cell in (u,v)-space + // or (x,y,z)-space. The point returned by GetCenterRaw is not necessarily + // unit length. + S2Point GetCenter() const; + S2Point GetCenterRaw() const; + + // Return the average area for cells at the given level. + static double AverageArea(int level); + + // Return the average area of cells at this level. This is accurate to + // within a factor of 1.7 (for S2_QUADRATIC_PROJECTION) and is extremely + // cheap to compute. + double AverageArea() const; + + // Return the approximate area of this cell. This method is accurate to + // within 3% percent for all cell sizes and accurate to within 0.1% for + // cells at level 5 or higher (i.e. squares 350km to a side or smaller + // on the Earth's surface). It is moderately cheap to compute. + double ApproxArea() const; + + // Return the area of this cell as accurately as possible. This method is + // more expensive but it is accurate to 6 digits of precision even for leaf + // cells (whose area is approximately 1e-18). + double ExactArea() const; + + // Return the distance from the cell to the given point. Returns zero if + // the point is inside the cell. + S1ChordAngle GetDistance(const S2Point& target) const; + + // Return the distance from the cell boundary to the given point. + S1ChordAngle GetBoundaryDistance(const S2Point& target) const; + + // Return the minimum distance from the cell to the given edge AB. Returns + // zero if the edge intersects the cell interior. + S1ChordAngle GetDistance(const S2Point& a, const S2Point& b) const; + + //////////////////////////////////////////////////////////////////////// + // S2Region interface (see s2region.h for details): + + S2Cell* Clone() const override; + S2Cap GetCapBound() const override; + S2LatLngRect GetRectBound() const override; + bool Contains(const S2Cell& cell) const override; + bool MayIntersect(const S2Cell& cell) const override; + + // Return true if the cell contains the given point "p". Note that unlike + // S2Loop/S2Polygon, S2Cells are considered to be closed sets. This means + // that points along an S2Cell edge (or at a vertex) belong to the adjacent + // cell(s) as well. + // + // If instead you want every point to be contained by exactly one S2Cell, + // you will need to convert the S2Cells to S2Loops (which implement point + // containment this way). + // + // The point "p" does not need to be normalized. + bool Contains(const S2Point& p) const override; +}; +``` + +See `s2cell.h` for additional methods. + +## S2CellUnion + +An `S2CellUnion` is an `S2Region` consisting of cells of various sizes. A cell +union is typically used to approximate some other shape. There is a tradeoff +between the accuracy of the approximation and how many cells are used. Unlike +polygons, cells have a fixed hierarchical structure. This makes them more +suitable for spatial indexing. + +An `S2CellUnion` is represented as a vector of sorted, non-overlapping +`S2CellIds`. By default the vector is also "normalized", meaning that groups +of 4 child cells have been replaced by their parent cell whenever possible. +S2CellUnions are not required to be normalized, but certain operations will +return different results if they are not (e.g., `Contains(S2CellUnion)`.) + +Here is the interface: + +```c++ +class S2CellUnion final : public S2Region { + public: + // Creates an empty cell union. + S2CellUnion() {} + + // Constructs a cell union with the given S2CellIds, then calls Normalize() + // to sort them, remove duplicates, and merge cells when possible. (See + // FromNormalized if your vector is already normalized.) + // + // The argument is passed by value, so if you are passing a named variable + // and have no further use for it, consider using std::move(). + // + // A cell union containing a single S2CellId may be constructed like this: + // + // S2CellUnion example({cell_id}); + explicit S2CellUnion(std::vector cell_ids); + + // Convenience constructor that accepts a vector of uint64. Note that + // unlike the constructor above, this one makes a copy of "cell_ids". + explicit S2CellUnion(const std::vector& cell_ids); + + // Constructs a cell union from S2CellIds that have already been normalized + // (typically because they were extracted from another S2CellUnion). + // + // The argument is passed by value, so if you are passing a named variable + // and have no further use for it, consider using std::move(). + // + // REQUIRES: "cell_ids" satisfies the requirements of IsNormalized(). + static S2CellUnion FromNormalized(std::vector cell_ids); + + // Constructs a cell union from a vector of sorted, non-overlapping + // S2CellIds. Unlike the other constructors, FromVerbatim does not require + // that groups of 4 child cells have been replaced by their parent cell. In + // other words, "cell_ids" must satisfy the requirements of IsValid() but + // not necessarily IsNormalized(). + // + // Note that if the cell union is not normalized, certain operations may + // return different results (e.g., Contains(S2CellUnion)). + // + // REQUIRES: "cell_ids" satisfies the requirements of IsValid(). + static S2CellUnion FromVerbatim(std::vector cell_ids); + + // Constructs a cell union that corresponds to a continuous range of cell + // ids. The output is a normalized collection of cell ids that covers the + // leaf cells between "min_id" and "max_id" inclusive. + // + // REQUIRES: min_id.is_leaf(), max_id.is_leaf(), min_id <= max_id. + static S2CellUnion FromMinMax(S2CellId min_id, S2CellId max_id); + + // Like FromMinMax() except that the union covers the range of leaf cells + // from "begin" (inclusive) to "end" (exclusive), as with Python ranges or + // STL iterator ranges. If (begin == end) the result is empty. + // + // REQUIRES: begin.is_leaf(), end.is_leaf(), begin <= end. + static S2CellUnion FromBeginEnd(S2CellId begin, S2CellId end); + + // Clears the contents of the cell union and minimizes memory usage. + void Clear(); + + // Gives ownership of the vector data to the client without copying, and + // clears the content of the cell union. The original data in cell_ids + // is lost if there was any. + std::vector Release(); + + // Convenience methods for accessing the individual cell ids. + int num_cells() const; + S2CellId cell_id(int i) const; + + // Vector-like methods for accessing the individual cell ids. + size_t size() const; + bool empty() const; + S2CellId operator[](int i) const; + + // Direct access to the underlying vector for STL algorithms. + const std::vector& cell_ids() const; + + // Returns true if the cell union is valid, meaning that the S2CellIds are + // valid, non-overlapping, and sorted in increasing order. + bool IsValid() const; + + // Returns true if the cell union is normalized, meaning that it is + // satisfies IsValid() and that no four cells have a common parent. + // Certain operations such as Contains(S2CellUnion) will return a different + // result if the cell union is not normalized. + bool IsNormalized() const; + + // Normalizes the cell union by discarding cells that are contained by other + // cells, replacing groups of 4 child cells by their parent cell whenever + // possible, and sorting all the cell ids in increasing order. Returns true + // if the number of cells was reduced. + bool Normalize(); + + // Replaces "output" with an expanded version of the cell union where any + // cells whose level is less than "min_level" or where (level - min_level) + // is not a multiple of "level_mod" are replaced by their children, until + // either both of these conditions are satisfied or the maximum level is + // reached. + // + // This method allows a covering generated by S2RegionCoverer using + // min_level() or level_mod() constraints to be stored as a normalized cell + // union (which allows various geometric computations to be done) and then + // converted back to the original list of cell ids that satisfies the + // desired constraints. + void Denormalize(int min_level, int level_mod, + std::vector* output) const; + + // If there are more than "excess" elements of the cell_ids() vector that + // are allocated but unused, reallocate the array to eliminate the excess + // space. This reduces memory usage when many cell unions need to be held + // in memory at once. + void Pack(int excess = 0); + + // Return true if the cell union contains the given cell id. Containment is + // defined with respect to regions, e.g. a cell contains its 4 children. + // This is a fast operation (logarithmic in the size of the cell union). + bool Contains(S2CellId id) const; + + // Return true if the cell union intersects the given cell id. + // This is a fast operation (logarithmic in the size of the cell union). + bool Intersects(S2CellId id) const; + + // Return true if this cell union contain/intersects the given other cell + // union. + bool Contains(const S2CellUnion& y) const; + bool Intersects(const S2CellUnion& y) const; + + // Returns the union of the two given cell unions. + S2CellUnion Union(const S2CellUnion& y) const; + + // Returns the intersection of the two given cell unions. + S2CellUnion Intersection(const S2CellUnion& y) const; + + // Specialized version of GetIntersection() that returns the intersection of + // a cell union with an S2CellId. This can be useful for splitting a cell + // union into pieces. + S2CellUnion Intersection(S2CellId id) const; + + // Returns the difference of the two given cell unions. + S2CellUnion Difference(const S2CellUnion& y) const; + + // Expands the cell union by adding a "rim" of cells on expand_level + // around the union boundary. + // + // For each cell c in the union, we add all cells at level + // expand_level that abut c. There are typically eight of those + // (four edge-abutting and four sharing a vertex). However, if c is + // finer than expand_level, we add all cells abutting + // c.parent(expand_level) as well as c.parent(expand_level) itself, + // as an expand_level cell rarely abuts a smaller cell. + // + // Note that the size of the output is exponential in + // "expand_level". For example, if expand_level == 20 and the input + // has a cell at level 10, there will be on the order of 4000 + // adjacent cells in the output. For most applications the + // Expand(min_radius, max_level_diff) method below is easier to use. + void Expand(int expand_level); + + // Expand the cell union such that it contains all points whose distance to + // the cell union is at most "min_radius", but do not use cells that are + // more than "max_level_diff" levels higher than the largest cell in the + // input. The second parameter controls the tradeoff between accuracy and + // output size when a large region is being expanded by a small amount + // (e.g. expanding Canada by 1km). For example, if max_level_diff == 4 the + // region will always be expanded by approximately 1/16 the width of its + // largest cell. Note that in the worst case, the number of cells in the + // output can be up to 4 * (1 + 2 ^ max_level_diff) times larger than the + // number of cells in the input. + void Expand(S1Angle min_radius, int max_level_diff); + + // The number of leaf cells covered by the union. + // This will be no more than 6*2^60 for the whole sphere. + uint64 LeafCellsCovered() const; + + // Approximate this cell union's area by summing the average area of + // each contained cell's average area, using the AverageArea method + // from the S2Cell class. + // + // NOTE: Since this is proportional to LeafCellsCovered(), it is + // always better to use the other function if all you care about is + // the relative average area between objects. + double AverageBasedArea() const; + + // Calculate this cell union's area by summing the approximate area for each + // contained cell, using the ApproxArea method from the S2Cell class. + double ApproxArea() const; + + // Calculate this cell union's area by summing the exact area for each + // contained cell, using the Exact method from the S2Cell class. + double ExactArea() const; + + // Return true if two cell unions are identical. + friend bool operator==(const S2CellUnion& x, const S2CellUnion& y); + + // Return true if two cell unions are different. + friend bool operator!=(const S2CellUnion& x, const S2CellUnion& y); + + //////////////////////////////////////////////////////////////////////// + // S2Region interface (see s2region.h for details): + + S2CellUnion* Clone() const override; + S2Cap GetCapBound() const override; + S2LatLngRect GetRectBound() const override; + bool Contains(const S2Cell& cell) const override; + bool MayIntersect(const S2Cell& cell) const override; + + // The point 'p' does not need to be normalized. + // This is a fast operation (logarithmic in the size of the cell union). + bool Contains(const S2Point& p) const override; + + // Appends a serialized representation of the S2CellUnion to "encoder". + void Encode(Encoder* const encoder) const; + + // Decodes an S2CellUnion encoded with Encode(). Returns true on success. + bool Decode(Decoder* const decoder); + +}; +``` + +## Approximating Regions + +An `S2RegionCoverer` is a class that allows arbitrary regions to be approximated +as unions of cells (`S2CellUnion`). This is useful for implementing various +sorts of search and precomputation operations. + +Typical usage: + + S2RegionCoverer::Options options; + options.set_max_cells(5); + S2RegionCoverer coverer(options); + S2Cap cap(center, radius); + S2CellUnion covering = coverer.GetCovering(cap); + +The result is a cell union of at most 5 cells that is guaranteed to cover the +given cap (a disc-shaped region on the sphere). + +The approximation algorithm is not optimal but does a pretty good job in +practice. The output does not always use the maximum number of cells allowed, +both because this would not always yield a better approximation, and because +`max_cells()` is a limit on how much work is done exploring the possible +covering as well as a limit on the final output size. + +The default `max_cells()` value of 8 is a reasonable point on the +precision/performance tradeoff curve for caps, but you may want to use higher +values for odd-shaped regions, such as long skinny lat-lng rects. [Here are some +examples of approximations.](examples/coverings) + +Here is the interface of `S2RegionCoverer`: + +```c++ +class S2RegionCoverer { + public: + class Options { + public: + // Sets the desired maximum number of cells in the approximation. Note + // that min_level() takes priority over max_cells(), i.e. cells below the + // given level will never be used even if this causes a large number of + // cells to be returned. + // + // Accuracy is measured by dividing the area of the covering by the area + // of the original region. The following table shows the median and worst + // case values for this area ratio on a test case consisting of 100,000 + // spherical caps of random size (generated using s2regioncoverer_test): + // + // max_cells: 3 4 5 6 8 12 20 100 1000 + // median ratio: 5.33 3.32 2.73 2.34 1.98 1.66 1.42 1.11 1.01 + // worst case: 215518 14.41 9.72 5.26 3.91 2.75 1.92 1.20 1.02 + // + // The default value of 8 gives a reasonable tradeoff between the number + // of cells used and the accuracy of the approximation. + // + // DEFAULT: kDefaultMaxCells + static constexpr int kDefaultMaxCells = 8; + int max_cells() const; + void set_max_cells(int max_cells); + + // Sets the minimum and maximum cell levels to be used. The default is to + // use all cell levels. + // + // To find the cell level corresponding to a given physical distance, use + // the S2Cell metrics defined in s2metrics.h. For example, to find the + // cell level that corresponds to an average edge length of 10km, use: + // + // int level = + // S2::kAvgEdge.GetClosestLevel(S2Earth::KmToRadians(length_km)); + // + // Note that min_level() takes priority over max_cells(), i.e. cells below + // the given level will never be used even if this causes a large number + // of cells to be returned. + // + // REQUIRES: max_level() >= min_level() + // DEFAULT: 0 + int min_level() const; + void set_min_level(int min_level); + + // DEFAULT: S2CellId::kMaxLevel + int max_level() const; + void set_max_level(int max_level); + + // Convenience function that sets both the maximum and minimum cell levels. + void set_fixed_level(int level); + + // If specified, then only cells where (level - min_level) is a multiple + // of "level_mod" will be used (default 1). This effectively allows the + // branching factor of the S2CellId hierarchy to be increased. Currently + // the only parameter values allowed are 1, 2, or 3, corresponding to + // branching factors of 4, 16, and 64 respectively. + // + // DEFAULT: 1 + int level_mod() const; + void set_level_mod(int level_mod); + }; + + // Constructs an S2RegionCoverer with the given options. + explicit S2RegionCoverer(const Options& options); + + // Default constructor. Options can be set using mutable_options(). + S2RegionCoverer(); + ~S2RegionCoverer(); + + // Returns the current options. Options can be modifed between calls. + const Options& options() const; + Options* mutable_options(); + + // Returns an S2CellUnion that covers (GetCovering) or is contained within + // (GetInteriorCovering) the given region and satisfies the current options. + // + // Note that if options().min_level() > 0 or options().level_mod() > 1, the + // by definition the S2CellUnion may not be normalized, i.e. there may be + // groups of four child cells that can be replaced by their parent cell. + S2CellUnion GetCovering(const S2Region& region); + S2CellUnion GetInteriorCovering(const S2Region& region); +}; +``` + +## Appendix: Alternatives Considered + +The S2 subdivision scheme has a number of advantages over the Hierarchical +Triangular Mesh (http://skyserver.org/HTM) framework: + +* Converting from a cell id to a point or vice versa is about 100 times faster + than the corresponding HTM operation. For example, a unit vector can be + converted to an S2CellId in about 0.15 microseconds, while converting a unit + vector to an HTM triangle id takes about 25 microseconds. Similarly, + converting an S2CellId to a vector takes about 0.04 microseconds, while the + same HTM operation takes about 19 microseconds. (If full resolution is not + needed, the HTM times can be improved by using fewer levels -- e.g. by using + 15 levels rather than 30, conversions are about twice as fast, but the + maximum resolution is reduced from 1cm to about 300 meters. The S2 library + is still about 100 times faster, though.) + +* The S2 library has no storage requirements. In contrast, the HTM library + precomputes the upper levels of the mesh during initialization to get even + the performance numbers mentioned above. Computing the first 6 levels (the + default) takes about 2MB of memory and 5 milliseconds. It is not practical + to precompute much more than this because memory requirements increase by a + factor of 4 for each additional level. + +* A linear scan of a set of S2CellIds follows a space-filling curve, which + maximizes locality of reference. For example, if each cell id is looked up + in a spatial index of some sort, a space-filling curve minimizes the number + of times that we switch from one index element to the next. This helps cache + performance (including non-local caches such as + [BigTable](https://cloud.google.com/bigtable){:target="_blank"}). A linear + scan through HTM triangles also has fairly good locality of reference due to + the hierarchical structure, but it does not define a space-filling curve + (due to the triangle ordering chosen by its inventors) and therefore the + locality of reference is not as good. (At any given level, the HTM path is + about 2.2 times longer than the corresonding S2CellId path.) + +Another alternative is HEALPix +[http://www.eso.org/science/healpix](http://www.eso.org/science/healpix){:target="_blank"}. +The main advantage of HEALPix is that it is suitable for calculations involving +spherical harmonics. This isn't relevant to any of our current applications, and +the scheme otherwise has several disadvantages from our point of view (e.g. cell +boundaries are not geodesics, base structure is more complicated). + +The COBE quadrilateralized spherical cube projection +(http://lambda.gsfc.nasa.gov/product/cobe/skymap_info_new.cfm) also starts by +projecting the six faces of a cube onto the unit sphere, and subdivides each +face in a quadtree fashion. However, it does not use a space-filling curve +scheme for labelling the cells, the cell edges are not geodesics, and it uses a +much more complicated projection scheme designed to minimize distortion. + +[^s2cell-locality]: Ideally the converse property would also be true (i.e., if + two cells are close together then their `S2CellIds` are also close + together), but unfortunately this turns out not to be possible. + +[^cell-center-parameters]: There are actually 3 different Hilbert curve + parameters for the point (0.5, 0.5), corresponding to the fact that the + Hilbert curve visits this point 3 times. The three parameters are + 0.1000..., 0.00101010..., and 0.11010101... (note that there is another + expansion 0.0111... which also visits this point, but it is numerically + equal to 0.1000...). This is similar to the fact that there are two + decimal expansions of the number 1 (1.000... and 0.999...). In any case, + the S2 library only uses the 0.1000... parameter for the cell center. diff --git a/docs/devguide/s2chain_interpolation_query.md b/docs/devguide/s2chain_interpolation_query.md new file mode 100644 index 00000000..9fce44d0 --- /dev/null +++ b/docs/devguide/s2chain_interpolation_query.md @@ -0,0 +1,228 @@ +--- +title: Interpolation along Chains +--- + +The purpose of this document is to provide a description of +[`S2ChainInterpolationQuery`](https://source.corp.google.com/piper///depot/google3/util/geometry/s2chain_interpolation_query.h) +and to illustrate its usage with code snippets. + +## Overview + +[`S2ChainInterpolationQuery`](https://source.corp.google.com/piper///depot/google3/util/geometry/s2chain_interpolation_query.h) +is a tool for interpolating point locations along chains. For example, it can +compute the position that is 50% along a polyline, or generate all the +positions on a polygon boundary that are 1e-4 radians apart. The interpolation +is based on spherical distance that is computed cumulatively along the edges of +the chain. It works with geometry classes that are derived from +[`S2Shape`](https://source.corp.google.com/piper///depot/google3/util/geometry/s2shape.h) +and represent polylines, polygons, or sets thereof. + +An `S2Shape` is organized as a collection of edges, which are identified by a +contiguous range of integer edge ids, starting at 0. The edges are grouped into +chains, where each chain consists of a sequence of edges connected end-to-end +and thus forming a polyline, or a closed polygon. For example, a shape +representing a single open polyline with _N_ vertices would consist of just one +chain containing _N - 1_ edges. A shape representing an island with a lake in it +would consist of two chains, one chain for the outer contour of the island and +the other for the contour of the lake. + +If the input shape contains multiple chains, the caller can specify the index of +the chain to be used, or let the query use all chains contained in the shape. +When multiple chains are used, they are assumed to be arranged end to end in the +distance space, in the same order as they are stored in the `S2Shape`, see +[Working with Multiple Chains](https://g3doc.corp.google.com/devguide/s2chain_interpolation_query.md#working-with-multiple-chains) +section for more details. + +`S2ChainInterpolationQuery` can be used for: + + * Finding the point that is located at a given distance from the first vertex + along the chain + * Finding the point that is located at a given normalized distance, that is a + fraction of the polyline's total length. + * Obtaining the edge ids of the resulting points + +## Initializing an `S2ChainInterpolationQuery` + +An instance of `S2ChainInterpolationQuery` can be initialized as follows: + +``` +const S2Shape* shape = ...; +int chain_id = ...; + +S2ChainInterpolationQuery query(shape, chain_id); +``` + +where `chain_id` is the non-negative id of the shape's chain on which the points +are going to be queried (0 <= `chain_id` < shape->num_chains()). + +If `chain_id` is negative (or altogether omitted) then all shape's chains are +going to be used in the order of their respective ids, as if there were a single +polyline. Note that in this case the resulting point may not be a continuous +function of distance, since the beginning of a chain may not coincide with the +end of the previous chain. + +When the input shape is known to contain only one chain, the `chain_id` +parameter can be either set to 0 or omitted. + +## Using the Query Methods + +Once an instance of `S2ChainInterpolationQuery` is initialized, the total length +of a chain (or chains) can be accessed via its `GetLength()` method: + +``` +S1Angle total_length = query.GetLength(); +``` + +Note that if the query is uninitialized, or initialized with an empty shape +(i.e. a shape that contains no chains), the returned total length is zero +(`S1Angle::Zero()`. + +To find the point at a given spherical distance _d_ along the chain(s), use the +`AtDistance()` method. This method also computes the id of the edge on which the +resulting point is located, as well as the actual parametric distance of the +point. The resulting point, edge id and distance can be obtained via the +respective accessor methods of `S2ChainInterpolationQuery::Result`: + +``` +S1Angle d = ...; +S2ChainInterpolationQuery::Result result = query.AtDistance(d); + +if (result.is_valid()) {} + S2Point point_at_distance = result.point(); + int edge_id = result.edge_id(); + S1Angle actual_distance = result.distance(); +} +``` + +If the input distance _d_ is less than or equal to the chain's total length, +then the actual distance returned by `result.distance()` equals _d_. If, +however, _d_ is greater than the chain's total length, then the resulting point +is snapped to the last vertex of the chain(s) and `result.distance()` returns +the total length of the chain. + +For example, consider the polyline _abc_ formed by three vertices _A_, _B_ and +_C_ with the respective latitudes/longitudes of (0, 0), (0, 1) and (0, 2). It +has two edges, _AB_ and _BC_, each of them having the spherical length of 1 +degree. For this polyline we have: + +``` +// Vertices A, B and C +auto vertices = s2textformat::ParsePointsOrDie("0:0, 0:1, 0:2"); + +S2LaxPolylineShape shape(vertices); + +// The above shape is known to contain just one chain (A,B,C), therefore one can +// omit the chain_id parameter in the query initialization. +S2ChainInterpolationQuery query(&shape); + +// The total length is 2 degrees = sum of lengths of edges (A,B) and (B,C). +S1Angle total_length = query.GetLength(); + +// Find the point at the spherical distance of 1.5 degrees. +S2ChainInterpolationQuery::Result result = query.AtDistance(S1Angle::Degrees(1.5))); +S2Point point = result.point(); +S1Angle distance = result.distance(); +int edge_id = result.edge_id(); +``` + +The resulting point is located half-way between vertices *B* and *C*, the +resulting distance is the same as the input distance of 1.5 degrees and the edge +id is 1, corresponding to the second edge of the chain. Alternatively, calling +`query.AtDistance(2.5)` yields the resulting point coinciding with *C*, a +resulting distance equal to the chain's total length (2 degrees), and an edge id +of 1, since the results of queries with distances exceeding the chain's length +are snapped to the last vertex of the chain. + +To find the point at a given normalized distance (i.e. at a fraction of the +total length) _f_, where 0 <= _f_ <= 1, use the `AtFraction()` method: + +``` +double f = ...; + +auto result = query.AtFraction(f); +S2Point point_at_distance = result.point(); +int edge_id = result.edge_id(); +S1Angle actual_distance = result.distance(); +``` + +The actual distance returned by `AtFraction()` is always within the interval of +[0, _total_length_]. If passed an input fraction value greater than 1, then +the resulting point is snapped to the last vertex of the chain(s) and the +distance returned by `AtFraction()` equals the total length of the chain. + +Calling `query.AtFraction(0.75)` for the polyline shape in the above +example would yield the resulting point half-way between _B_ and _C_, resulting +distance of 1.5 and the edge id = 1 (the last edge of the chain). And +`query.AtFraction(1.5)` gives the resulting point coinciding with the +last vertex _C_, resulting length = 2 (the total length) and the resulting edge +id = 1 (the edge corresponding to the last vertex). + +Both `AtDistance()` and `AtFraction()` return valid results for all input +distance values if and only if the query has been initialized with a shape +containing at least one valid edge. + +Let us now consider an +[`S2LaxPolygonShape`](https://source.corp.google.com/piper///depot/google3/util/geometry/s2lax_polygon_shape.h) +containing two chains, each being 1 degree long and containing two edges. + +If the query for this shape is initialized with a `chain_id` of 0: + +``` +S2ChainInterpolationQuery query(&shape, 0); +``` + +then only the first chain (chain_id = 0) is considered and the other chain is +ignored. In this case `query.GetLength()` returns the value corresponding to the +angle of 1 degree, which is the length of the single chain. + +Now, if the query is initialized with a `chain_id` that is negative (or +omitted), both chains are considered as if they were part of the same polyline, +and `query.GetLength()` returns an angle of 2 degrees, which is the length of +both chains combined: + +``` +S2ChainInterpolationQuery query(&shape); +S1Angle total_length = query.GetLength(); // total_length.degrees() == 2 +``` + +In this case calling `query.AtDistance(1.5)` or `query.AtFraction(0.75)` would +yield a point halfway down the second chain, and the resulting an edge id = 1 +(the second edge of the shape corresponding to the first and the only edge of +the second chain). + +## Working with Multiple Chains + +An important caveat when considering multiple chain is that the resulting point +is not a continuous function of the input distance, since the last vertex of a +chain may not coincide with the first vertex of the next chain. In the above +case, calling `query.AtFraction( 0.5 )` yields the resulting point = last vertex +of the first chain, while `query.AtFraction( 0.5 + epsilon )` yields the first +vertex of the second chain, which may be located far apart from the last vertex +of the previous chain. + +A typical use case for using queries with multiple chains is when one needs to +generate an *approximately* uniform sampling of all chains contained in a shape. +In these cases one can just call `query.AtFraction(t)` with *t* ranging from 0.0 +to 1.0 with a certain step, depending on the desired sampling density. Note that +such sampling is only approximately uniform, since the chain lengths won't +necessarily be multiples of the step size. + +## Computational Complexity + +* `S2ChainInterpolationQuery` initialization has the computational complexity + of *O(N)*, where *N* is the total number of edges. +* `AtDistance()` and `AtFraction()` calls have the complexity of *O( log(N) + )*. +* `GetLength()` call has the constant complexity *O(1)*. +* The memory footprint of the query is *O(N)*. + +## Numerical Accuracy + +The `S2Point` returned by `result.point()` accessor is an approximation of the +true position with the error bound by +[`kGetPointOnLineError`](https://source.corp.google.com/piper///depot/google3/util/geometry/s2edge_distances.h;l=148?q=kGetPointOnLineError). + +The computation error of the total length is bound by 3.25 * DBL_EPSILON * _N_, +where _N_ is the number of edges. An additional loss of precision may occur due +to floating point multiplication of `fraction * total_length` during +`AtFraction()` call. diff --git a/docs/devguide/s2closestedgequery.md b/docs/devguide/s2closestedgequery.md new file mode 100644 index 00000000..a3a21780 --- /dev/null +++ b/docs/devguide/s2closestedgequery.md @@ -0,0 +1,373 @@ +--- +title: Finding Nearby Edges +--- + +## Overview + +`S2ClosestEdgeQuery` is a tool for finding the closest edge(s) between two +geometries. By using the appropriate options, this class can answer questions +such as: + + - Find the minimum distance between two geometries *A* and *B*. + + - Find all edges of geometry *A* that are within a distance *d* of geometry + *B*. + + - Find the *k* edges of geometry *A* that are closest to a given point *P*. + +The input geometries may consist of any number of points, polylines, +and polygons (collectively known as *shapes*). Shapes do not need to +be disjoint; they may overlap or intersect arbitrarily. The +implementation is designed to be fast for both simple and complex +geometries (e.g., one edge vs. 100 million edges). + +## The Index + +`S2ClosestEdgeQuery` operates on two geometries, known as the *index* +and the *target*. Each query returns the subset of edges in the index +that have a specified relationship to the target (e.g., the closest +edge to the target, or all edges within a distance *d* of the target). + +The indexed geometry is provided as an `S2ShapeIndex` object, which is +simply a collection of points, polylines, and/or polygons (possibly +intersecting, as mentioned above). Here is a simple example showing +how to build an index containing several polylines: + + vector polylines = ...; + S2ShapeIndex index; + for (S2Polyline* polyline : polylines) { + index.Add(new S2Polyline::Shape(polyline)); + } + +An `S2ClosestEdgeQuery` object can then be constructed like so: + + S2ClosestEdgeQuery query(&index); + +Note that `index` is passed as a const pointer rather than a reference +to reflect the fact that the index must persist for the lifetime of the +query object. + +## The Target + +The *target* represents the geometry that distances are measured to. +It is specified as an object of type `S2ClosestEdgeQuery::Target`. +There are `Target` subtypes for measuring the distance to a point, an +edge, an `S2Cell`, or an entire `S2ShapeIndex` (which allows finding +closest edges or measuring distances between two arbitrary collections +of geometry). + +For example, the following code tests whether any of the polylines in +the example above is closer than `dist_limit` to a given point +`test_point`: + + S2Point test_point = ...; + S1ChordAngle dist_limit = ...; + S2ClosestEdgeQuery::PointTarget target(test_point); + if (query.IsDistanceLess(&target, dist_limit)) { ... } + +Note that `target` is passed as a pointer; this is because `Target` +objects may contain temporary data and can therefore only be used by +one query at a time. (The target and query object are intended for use +only within a single thread, as opposed to the `S2ShapeIndex` object +which can be used in many different threads simultaneously.) + +## Query Methods + +The fastest and simplest query method is `IsDistanceLess`, which +returns true if the distance to the target is less than a given limit: + + bool IsDistanceLess(Target* target, S1ChordAngle limit); + +If you want to test whether the distance is less than or equal to "limit" +instead, you can use `IsDistanceLessOrEqual` instead: + + if (query.IsDistanceLessOrEqual(&target, limit)) { ... } + +To compute the actual distance, you can use the `GetDistance` method: + + S1ChordAngle GetDistance(Target* target); + +However note that if you only need to compare the distance against some +threshold value, it is much faster to call `IsDistanceLess` instead. + +The return type of this method (`S1ChordAngle`) is an efficient +representation of the distance between two points on the unit sphere. +It can easily be converted to an `S1Angle` if desired (which measures +angles in radians or degrees), and from there to an actual distance on +the Earth spheroid. (See also [Modeling Accuracy](#modeling-accuracy).) + +The most general query method is `FindClosestEdges`: + + std::vector FindClosestEdges(Target* target); + +This method finds the closest edge(s) to the target that satisfy the +query options (described in the following section). Each edge is +returned in the form of a `Result` object that identifies the edge and +specifies the distance to that edge: + + class Result { + Distance distance() const; // The distance from the target to this edge. + int32 shape_id() const; // Identifies an indexed shape. + int32 edge_id() const; // Identifies an edge within the shape. + }; + +The `Distance` type is simply a thin wrapper around `S1ChordAngle` that +provides some additional methods needed by the `S2ClosestEdgeQuery` +implementation. The other two fields identify an edge in the +`S2ShapeIndex` that was passed to the `S2ClosestEdgeQuery` +constructor. The query object has two convenience methods that may be +useful for processing result edges: + + // Returns the endpoints of the given result edge. + S2Shape::Edge GetEdge(const Result& result) const; + + // Returns the point on given result edge that is closest to "point". + S2Point Project(const S2Point& point, const Result& result) const; + +For example, to retrieve the endpoints of each edge and also find the +closest point on each result edge to a given target point, you could do +this: + + for (const auto& result : query.FindClosestEdges(&target)) { + S2Shape::Edge edge = query.GetEdge(result); + S2Point closest_point = query.Project(point, result); + } + +There is also a convenience method for the common case where only one +edge (the closest edge is needed). This method returns a single +`Result` object rather than a vector: + + Result FindClosestEdge(Target* target); + +## Options + +`S2ClosestEdgeQuery` supports a number of options that control which +edges are returned by each query. These options can be specified as a +second argument to the constructor: + + S2ClosestEdgeQuery::Options options; + S2ClosestEdgeQuery query(&index, options); + +Unlike most classes in the library, `S2ClosestEdgeQuery` also allows +options to be modified between calls to query methods via the +`mutable_options()` accessor method. For example: + + S1ChordAngle distance1 = query.GetDistance(&target1); + query.mutable_options()->set_max_distance(dist_limit); + S1ChordAngle distance2 = query.GetDistance(&target2); + +This can be significantly more efficient compared to changing the +options by constructing a new `S2ClosestEdgeQuery` object, because the +query object contains a significant amount of internal state that is +reused from one query to the next. + +The most important options are the following: + + // Specifies that at most "max_results" edges should be returned. + // + // REQUIRES: max_results >= 1 + // DEFAULT: numeric_limits::max() + int max_results() const; + void set_max_results(int max_results); + + // Specifies that only edges whose distance to the target is less than + // "max_distance" should be returned. + // + // Note that edges whose distance is exactly equal to "max_distance" are + // not returned. Normally this doesn't matter, because distances are not + // computed exactly in the first place, but if such edges are needed then + // see set_inclusive_max_distance() below. + // + // DEFAULT: Distance::Infinity() + void set_max_distance(S1ChordAngle max_distance); + + // Like set_max_distance(), except that edges whose distance is exactly + // equal to "max_distance" are also returned. Equivalent to calling + // set_max_distance(max_distance.Successor()). + void set_inclusive_max_distance(S1ChordAngle max_distance); + +Note that you will always want to set one of these options before +calling `FindClosestEdges`, because otherwise all the edges in the +entire `S2ShapeIndex` will be returned. (This is not necessary +when calling `FindClosestEdge`, `GetDistance`, or `IsDistanceLess`, +because these methods all implicitly limit max_results() to 1.) + +Another important option is `include_interiors()`, which determines +whether distances are measured to the boundary and interior of polygons +or only to their boundaries: + + // Specifies that polygon interiors should be included when measuring + // distances. + // + // DEFAULT: true + bool include_interiors() const; + void set_include_interiors(bool include_interiors); + +By default this option is true, which means for example that a point +inside a polygon has a distance of zero. Note that in this situation +there is no "closest edge" to return (since the minimum distance is +attained at a point interior to the polygon), and therefore the +`FindClosestEdge(s)` methods indicate this in the `Result` object(s) by +setting `edge_id` to -1. (Such results can be identified using the +`is_interior()` method.) For example: + + for (const auto& result : query.FindClosestEdges(&target)) { + if (result.is_interior()) { + ProcessInterior(result.shape_id()); + } else { + ProcessEdge(result.shape_id, result.edge_id); + } + } + +The following option can be useful when performance is critical: + + // Specifies that edges up to max_error() further away than the true + // closest edges may be substituted in the result set, as long as such + // edges satisfy all the remaining search criteria (such as max_distance). + // + // DEFAULT: Distance::Zero() + Distance max_error() const; + void set_max_error(Distance max_error); + +This option give the algorithm permission to stop the search early when +the best possible improvement in distance drops below `max_error()`. +For example, this option is used internally by the `IsDistanceLess` +predicate to stop the search as soon as it finds any edge whose +distance is less than the given limit, rather than continuing to search +for an edge that is even closer (which would be wasted effort). + +Note that this options only has an effect if `max_results()` is also +specified; otherwise all edges closer than `max_distance()` will always +be returned. + +## Target Types + +The following target types exist in addition to the `PointTarget` type +described earlier. + +`EdgeTarget` computes the closest distance to a spherical geodesic edge. For +example, to find the closest indexed edge to a target edge `ab`, you could do +this: + + S2ClosestEdgeQuery::EdgeTarget target(a, b); + auto result = query.FindClosestEdge(&target); + +`CellTarget` computes the distance to an `S2Cell`, including the interior of the +cell. This target type is used extensively internally to implement other +algorithms. For example, to test whether the distance to `cell` is less than +`limit`: + + S2ClosestEdgeQuery::CellTarget target(cell); + if (query.IsDistanceLess(&target, limit)) { ... } + +Finally, `ShapeIndexTarget` computes the distance to an `S2ShapeIndex` +(containing an arbitrary collection of points, polylines, and polygons, possibly +overlapping). For example, to calculate the distance between two +`S2ShapeIndexes`, you can do this: + + S2ClosestEdgeQuery query(&index1); + S2ClosestEdgeQuery::ShapeIndexTarget target(&index2); + S1ChordAngle distance = query.GetDistance(&target); + +(As always, if you only want to compare the distance against a threshold value +then it is much more efficient to call `IsDistanceLess` instead.) + +The `ShapeIndexTarget` type also supports its own `Options`. In particular, if +you want to measure the distance to the boundaries of polygons in the target +index (ignoring the polygon interiors), you can call + + target.set_include_interiors(false); + +before executing the distance query. Note that the `include_interiors()` option +on the query and the target are independent, so for example you can measure +distance from the boundary of one polygon to the boundary and interior of +another. + +## Modeling Accuracy + +This class models the Earth as a sphere, which can lead to distance errors of up +to 0.56% compared with modeling the Earth as an ellipsoid. While this is +perfectly adequate for many purposes (e.g., determining nearby restaurants), +some applications need higher accuracy. + +This can be achieved by recomputing distances using using an external geodesy +library such as `geographiclib`. For exact results, all distances would need to +be computed this way (and in fact there is a templatized base class +`S2ClosestEdgeQueryBase` that allows this), but in most cases it is sufficient +to use spherical distances to determine which pair of points is closest, and +then recalculate the distance accurately for that pair of points only. (This is +because the spherical distance error varies relatively slowly over the Earth's +surface.) For example, if the target is a point then you could do the +following: + + S2ClosestEdgeQuery::PointTarget target(target_point); + auto result = query.FindClosestEdge(&target); + double meters; + if (result.is_empty()) { + // The index is empty. + meters = std::numeric_limits::infinity(); + } else if (result.is_interior()) { + // The target point is inside an indexed shape. + meters = 0.0; + } else { + // Find the point on the result edge that is closest to the target. + S2Point closest_point = query.Project(target_point, result); + meters = GeoidDistance(target_point, closest_point); + } + +The same technique can be used for more complex targets (`ShapeIndexTarget`) by +executing a second query to find the closest pair of edges between the two +geometry collections, finding the points on that edge pair that achieve the +minimum spherical distance (see `S2::GetEdgePairClosestPoints`), and then +calculating the ellipsoidal distance between those points: + + S2ClosestEdgeQuery query1(&index1); + S2ClosestEdgeQuery::ShapeIndexTarget target2(&index2); + auto result1 = query1.FindClosestEdge(&target2); + double meters; + ... like the code above ... + } else { + // Get the edge from index1 (edge1) that is closest to index2. + S2Shape::Edge edge1 = query1.GetEdge(result1); + + // Now find the edge from index2 (edge2) that is closest to edge1. + S2ClosestEdgeQuery query2(&index2); + S2ClosestEdgeQuery::EdgeTarget target1(edge1.v0, edge1.v1); + auto result2 = query2.FindClosestEdge(&target1); + S2Shape::Edge edge2 = query2.GetEdge(result2); + + // Find the closest point pair on edge1 and edge2. + auto closest = S2::GetEdgePairClosestPoints(edge1.v0, edge1.v1, + edge2.v0, edge2.v1); + meters = GeoidDistance(closest.first, closest.second); + } + +## Numerical Accuracy + +The spherical geodesic distance calculation used by this class is fast and +accurate, but not exact. Expressing the errors in terms of distances on the +Earth's surface, for distances up to 10,000 km the error is less than 6 +nanometers. However for points that are nearly antipodal, the error can be as +large as 50 cm (dropping off to less than 1 micrometer for points that are 50 km +away from being antipodal). + +If you would like to measure distances conservatively by including these error +bounds, you can use the following option: + + options.set_conservative_max_distance(limit); + +This will automatically increase `limit` by the maximum error in the distance +calculation, ensuring that all edges whose true, exact distance is less than or +equal to `limit` are returned (along with some edges whose true distance is +slightly greater). + +Algorithms that need to do exact distance comparisons can use this option to +compute a set of candidate edges that can then be filtered further using exact +predicates (such as `s2pred::CompareEdgeDistance`). + +There is also a conservative version of the `IsDistanceLessOrEqual` predicate, + + bool IsConservativeDistanceLessOrEqual(Target* target, S1ChordAngle limit); + +which automatically increases the given limit by the maximum error in the +distance calculation. diff --git a/docs/devguide/s2shape_nesting_query.md b/docs/devguide/s2shape_nesting_query.md new file mode 100644 index 00000000..62e8fe1d --- /dev/null +++ b/docs/devguide/s2shape_nesting_query.md @@ -0,0 +1,99 @@ +--- +title: Classifying shape chains +--- + +## Overview + +`S2ShapeNestingQuery` is a query to divide the chains of an `S2Shape` into +either "shells", which are analogous to exterior rings in non-spherical +geometry, or "holes" which are contained by a shell. + +Shells never have a parent of their own, and may have zero or more holes in +them, while holes always have a single parent and never have holes of their own. + +Since we operate in spherical coordinates, the relationship between chains is +ambiguous. If you imagine a band around the equator, we can consider the +northern chain to be a shell and the southern chain to be a hole: + +![Band around the equator with northern edge marked as the shell](./img/nesting-equator-rings0.png) + +Or, equivalently, we can reverse that ordering and consider the southern edge to +be the shell: + +![Band around the equator with southern edge marked as the shell](./img/nesting-equator-rings1.png) + +To resolve this ambiguity, we have to specify one chain of the shape to be a +shell by definition, with the other chains classified in relation to it. We +refer to this shell as the *datum* shell, and the exact method we use to choose +it is given by an `S2DatumStrategy` which can be provided through the +`S2ShapeNestingQueryOptions` passed to the query. The default strategy is to use +the *first* chain in the shape as the datum shell. + +## Usage + +Assuming we have an `S2ShapeIndex` called `index`, and we're interested in the +shape at index 3: + +```c++ +S2ShapeNestingQuery nesting_query(&index); +vector relations = nesting_query.ShapeNesting(3); +``` + +The entries of `relations` are a 1:1 mapping with the chains in the shape. Thus, +`relations[2]` corresponds to chain 2 in shape 3. The `ChainRelation` class +itself provides information on the relationships between a given chain and other +chains, and whether a chain is a shell or a hole. + +As a concrete example, GeoJSON requires that each polygon have one shell and +zero or more holes. With the chain relationship information, we can easily scan +through our geometry and write it out as GeoJSON properly: + +```c++ +S2ShapeNestingQuery nesting_query(&index); +vector relations = nesting_query.ShapeNesting(3); + +for (int chain = 0; chain < relations.size(); chain++) + const ChainRelation relation& = relations[chain]; + if (relation.IsShell()) { + auto polygon = GeoJsonPolygon(chain) + for (int hole : relation.Holes()) { + polygon.AddHole(hole); + } + WriteGeoJsonPolygon(polygon); + } +} +``` + + +To perform chain classification, we use the fact that we have one chain that we +know *a priori* is a shell (via the datum strategy). We take a point on that +shell (the exact one doesn't matter), and draw a line to a random point on each +of the other chains. + +We look up the edges that are crossed by that line and keep track of potential +parents as a bit vector with one entry per chain. Every time we cross an edge, +we determine the chain it belongs to, and toggle that bit. + +When we've done that, we may have multiple possible parents for a given chain, +so we have to do a reduction step which looks like this: + +``` +for each chain0 with exactly one parent: + for each chain1: + if both chain0,parent in parents of chain1: + unset parent in parents of chain1 +``` + +This says that, if we know that A ⊃ B, and B ⊃ C, then C shouldn't *directly* +consider A as a possible parent. + +Lastly, we determine each chain's depth from its root parent, and chains that +are at an even depth have their parent link broken and thus become shells. + +To help visualize, this animation shows a test running to classify four very +concave nested chains. Chain 2 (blue) is being classified by drawing edges from +chain 0 (purple). The edge between the two points is show in red and intersected +edges are highlighted. To test the robustness of the algorithm, the test is run +multiple times, rotating the order of the vertices in the chains each time. + +![Animation of nesting query algorithm](./img/nesting-animation.gif) diff --git a/docs/devguide/s2shapeindex.md b/docs/devguide/s2shapeindex.md new file mode 100644 index 00000000..c2553bab --- /dev/null +++ b/docs/devguide/s2shapeindex.md @@ -0,0 +1,685 @@ +--- +title: The S2ShapeIndex Class +--- + +## Overview + +`S2ShapeIndex` is an abstract base class for indexing polygonal geometry +in memory. The objects in the index are known as "shapes", and may +consist of points, polylines, and/or polygons, possibly overlapping. The +index makes it very fast to answer queries such as finding nearby shapes, +measuring distances, testing for intersection and containment, etc. + +Each object in the index implements the `S2Shape` interface. An `S2Shape` +is a collection of edges that optionally defines an interior. The edges do +not need to be connected, so for example an `S2Shape` can represent a +polygon with multiple shells and/or holes, or a set of polylines, or a set +of points. All geometry within a single `S2Shape` must have the same +dimension, so for example if you want to create an `S2ShapeIndex` containing +a polyline and 10 points, then you will need at least two different shapes. + +There are a number of built-in classes that work with `S2ShapeIndex` +objects. Generally these classes are designed to work with any collection +of geometry that can be represented by an `S2ShapeIndex`, i.e. any +combination of points, polylines, and polygons. Such classes include: + +* [`S2ContainsPointQuery`](#s2containspointquery): returns the shape(s) + that contain a given point. +* [`S2ClosestEdgeQuery`](s2closestedgequery): returns the closest edge(s) + to a given point, edge, `S2CellId`, or `S2ShapeIndex`. +* `S2CrossingEdgeQuery`: returns the edge(s) that cross a given edge. +* `S2BooleanOperation`: computes boolean operations such as union, and + boolean predicates such as containment. +* `S2ShapeIndexRegion`: computes approximations for a collection of + geometry. +* `S2ShapeIndexBufferedRegion`: computes approximations that have been + expanded by a given radius. + +## Example + +Here is an example showing how to index a set of polygons and then +determine which polygon(s) contain each of a set of query points: + + void TestContainment(const vector& points, + const vector& polygons) { + MutableS2ShapeIndex index; + for (auto polygon : polygons) { + index.Add(absl::make_unique(polygon)); + } + auto query = MakeS2ContainsPointQuery(&index); + for (const auto& point: points) { + for (S2Shape* shape : query.GetContainingShapes(point)) { + S2Polygon* polygon = polygons[shape->id()]; + ... do something with (point, polygon) ... + } + } + } + +This example uses `S2Polygon::Shape`, which is one example of an `S2Shape` +object. `S2Polyline` and `S2Loop` also have nested `Shape` classes, and +there are additional `S2Shape` types defined in `*_shape.h`. + +## S2Shape + +The purpose of `S2Shape` is to represent polygonal geometry in a flexible +way. It is organized as a collection of edges that optionally defines an +interior. All geometry represented by an `S2Shape` must have the same +dimension, which means that an `S2Shape` can represent either a set of +points, a set of polylines, or a set of polygons. + +`S2Shape` is defined as an abstract base class in order to give clients +control over the underlying data representation. Sometimes an `S2Shape` +does not have any data of its own, but instead "wraps" some other class. +For example, `S2Polygon::Shape` wraps an `S2Polygon` object. Here are a +few commonly used types of `S2Shape`: + +`S2Polygon::Shape` +: a polygon + +`S2LaxPolygonShape` +: a polygon, possibly with degeneracies + +`S2Polyline::Shape` +: a polyline + +`S2LaxPolylineShape` +: a polyline, possibly with degeneracies + +`S2PointVectorShape` +: a collection of points + +`S2EdgeVectorShape` +: a collection of edges + +Some classes such as `S2Polygon` have a nested `Shape` class, and others +are defined as standalone objects. It is also easy for clients to +implement their own subtypes, since the `S2Shape` interface is fairly +minimal (see below). + +The edges of an `S2Shape` are identified by a contiguous range of *edge +ids* starting at 0. The edges are further subdivided into *chains*, where +each chain consists of a sequence of edges connected end-to-end (a +polyline). For example, an `S2Shape` representing two polylines AB and CDE +would have three edges (AB, CD, DE) grouped into two chains: (AB) and (CD, +DE). Similarly, an `S2Shape` representing 5 points would have 5 chains +consisting of one edge each. + +Each chain is represented by a contiguous range of edge ids with no gaps +between chains. Edges can be accessed using either the global numbering +(edge id) or a (chain_id, offset) pair. The global numbering is sufficient +for most purposes, but the chain representation is useful for some +algorithms such as polygon operations (see `S2BooleanOperation`). + +Here are the most important methods of the `S2Shape` interface: + +```c++ +class S2Shape { + public: + // An edge, consisting of two vertices "v0" and "v1". Zero-length edges are + // allowed, and can be used to represent points. + struct Edge { + S2Point v0, v1; + Edge(); + Edge(const S2Point& _v0, const S2Point& _v1); + friend bool operator==(const Edge& x, const Edge& y); + friend bool operator<(const Edge& x, const Edge& y); + }; + + S2Shape(); + virtual ~S2Shape(); + + // Returns the number of edges in this shape. Edges have ids ranging from 0 + // to num_edges() - 1. + virtual int num_edges() const = 0; + + // Returns the endpoints of the given edge id. + // + // REQUIRES: 0 <= id < num_edges() + virtual Edge edge(int edge_id) const = 0; + + // Returns the dimension of the geometry represented by this shape. + // + // 0 - Point geometry. Each point is represented as a degenerate edge. + // + // 1 - Polyline geometry. Polyline edges may be degenerate. A shape may + // represent any number of polylines. Polylines edges may intersect. + // + // 2 - Polygon geometry. Edges should be oriented such that the polygon + // interior is always on the left. In theory the edges may be returned + // in any order, but typically the edges are organized as a collection + // of edge chains where each chain represents one polygon loop. + // Polygons may have degeneracies (e.g., degenerate edges or sibling + // pairs consisting of an edge and its corresponding reversed edge). + // + // Note that this method allows degenerate geometry of different dimensions + // to be distinguished, e.g. it allows a point to be distinguished from a + // polyline or polygon that has been simplified to a single point. + virtual int dimension() const = 0; + + // A unique id assigned to this shape by S2ShapeIndex. Shape ids are + // assigned sequentially starting from 0 in the order shapes are added. + int id() const { return id_; } +}; +``` + +The additional methods related to edge chains are summarized below: + +```c++ + // A range of edge ids corresponding to a chain of zero or more connected + // edges, specified as a (start, length) pair. The chain is defined to + // consist of edge ids {start, start + 1, ..., start + length - 1}. + struct Chain { + int32 start, length; + Chain(); + Chain(int32 _start, int32 _length); + }; + + // The position of an edge within a given edge chain, specified as a + // (chain_id, offset) pair. Chains are numbered sequentially starting from + // zero, and offsets are measured from the start of each chain. + struct ChainPosition { + int32 chain_id, offset; + ChainPosition(); + ChainPosition(int32 _chain_id, int32 _offset); + }; + + // Returns the number of contiguous edge chains in the shape. For example, + // a shape whose edges are [AB, BC, CD, AE, EF] would consist of two chains + // (AB,BC,CD and AE,EF). Every chain is assigned a "chain id" numbered + // sequentially starting from zero. + // + // Note that it is always acceptable to implement this method by returning + // num_edges() (i.e. every chain consists of a single edge), but this may + // reduce the efficiency of some algorithms. + virtual int num_chains() const = 0; + + // Returns the range of edge ids corresponding to the given edge chain. The + // edge chains must form contiguous, non-overlapping ranges that cover the + // entire range of edge ids. This is spelled out more formally below: + // + // REQUIRES: 0 <= i < num_chains() + // REQUIRES: chain(i).length >= 0, for all i + // REQUIRES: chain(0).start == 0 + // REQUIRES: chain(i).start + chain(i).length == chain(i+1).start, + // for i < num_chains() - 1 + // REQUIRES: chain(i).start + chain(i).length == num_edges(), + // for i == num_chains() - 1 + virtual Chain chain(int chain_id) const = 0; + + // Returns the edge at offset "offset" within edge chain "chain_id". + // Equivalent to "shape.edge(shape.chain(chain_id).start + offset)" + // but may be more efficient. + virtual Edge chain_edge(int chain_id, int offset) const = 0; + + // Finds the chain containing the given edge, and returns the position of + // that edge as a (chain_id, offset) pair. + // + // REQUIRES: shape.chain(pos.chain_id).start + pos.offset == edge_id + // REQUIRES: shape.chain(pos.chain_id + 1).start > edge_id + // + // where pos == shape.chain_position(edge_id). + virtual ChainPosition chain_position(int edge_id) const = 0; +``` + +### S2Shape vs. S2Region + +It is worth pointing out that S2 defines two extensible interfaces for +representing geometry: `S2Shape` and `S2Region`. They can be compared as +follows: + +* The purpose of `S2Shape` is to flexibly represent polygonal geometry. + (This includes not only polygons, but also points and polylines.) Most of + the core S2 operations will work with any class that implements the + `S2Shape` interface. + +* The purpose of `S2Region` is to compute approximations for geometry. For + example, there are methods for computing bounding rectangles and discs, + and `S2RegionCoverer` can be used to approximate a region to any desired + accuracy as a collection of cells. Unlike `S2Shape`, an `S2Region` may + represent non-polygonal geometry such as discs (`S2Cap`). + +## MutableS2ShapeIndex + +The most important subtype of `S2ShapeIndex` is `MutableS2ShapeIndex`. +This class allows not only building an index but also updating it +incrementally by adding or removing shapes. Its major methods are as +follows: + +```c++ +class MutableS2ShapeIndex final : public S2ShapeIndex { + public: + // Create an S2ShapeIndex that uses the default option settings. + S2ShapeIndex(); + ~S2ShapeIndex(); + + // Takes ownership of the given shape and adds it to the index. Also assigns + // a unique id to the shape (shape->id()) and returns that id. Shape ids + // are assigned sequentially starting from 0 in the order shapes are added. + // Invalidates all iterators and their associated data. + int Add(std::unique_ptr shape); + + // Removes the given shape from the index and return ownership to the caller. + // Invalidates all iterators and their associated data. + std::unique_ptr Release(int shape_id); + + // Resets the index to its original state and returns ownership of all + // shapes to the caller. This method is more efficient than removing + // all shapes one at a time. + std::vector> ReleaseAll(); + + // Resets the index to its original state and deletes all shapes. + void Clear(); + + // The number of distinct shape ids that have been assigned. This equals + // the number of shapes in the index provided that no shapes have ever been + // removed. (Shape ids are not reused.) + int num_shape_ids() const override; + + // Returns a pointer to the shape with the given id, or nullptr if the shape + // has been removed from the index. + S2Shape* shape(int id) const override; + + // Returns the number of bytes currently occupied by the index (including any + // unused space at the end of vectors, etc). It has the same thread safety + // as the other "const" methods (see introduction). + size_t SpaceUsed() const; + + // Minimizes memory usage by requesting that any data structures that can be + // rebuilt should be discarded. This method invalidates all iterators. + // + // Like all non-const methods, this method is not thread-safe. + void Minimize() override; +}; +``` + +One important feature to be aware of is that updates are batched together +and applied lazily on the first subsequent query. This is important for +efficiency, and also means that we can avoid building the index entirely +when in fact there is no subsequent query. + +Note that even though updates are applied lazily, `MutableS2ShapeIndex` has +the same thread-safety properties as `vector`: const methods are +thread-safe, while non-const methods are not thread-safe. In other words, +no additional locking is required for multi-threaded read access. + +## S2ContainsPointQuery + +`S2ContainsPointQuery` is an example of a class that is designed to +work with `S2ShapeIndex`. It has methods to find the shapes and/or +edges that contain a given point. Here is an example showing how to +use it: + + S2ContainsPointQueryOptions options(S2VertexModel::CLOSED); + auto query = MakeS2ContainsPointQuery(&index, options); + if (query.Contains(point)) { + ... at least one shape contains "point" ... + } + +Other methods of `S2ContainsPointQuery` allow finding the set of shapes +containing a given point, visiting the edges incident to a given point, +etc. + +Note that because most of its methods are very fast, `S2ContainsPointQuery` +takes the `S2ShapeIndex` type as a template argument. This allows the +overhead of virtual method calls to be eliminated when the actual subtype +is known. `MakeS2ContainsPointQuery` is a convenience function that +constructs an instance appropriate for the given `S2ShapeIndex` argument: + +```c++ +template +inline S2ContainsPointQuery MakeS2ContainsPointQuery( + const IndexType* index, + const S2ContainsPointQueryOptions& options = S2ContainsPointQueryOptions()); +``` + +By default, polygon boundaries are modeled as *semi-open*, which means that +if several polygons tile the region around a shared vertex, then exactly +one of those polygons contains that vertex. If other behavior is desired, +it can be specified as follows: + +```c++ +// Defines whether shapes are considered to contain their vertices. Note that +// these definitions differ from the ones used by S2BooleanOperation. +// +// - In the OPEN model, no shapes contain their vertices (not even points). +// Therefore Contains(S2Point) returns true if and only if the point is +// in the interior of some polygon. +// +// - In the SEMI_OPEN model, polygon point containment is defined such that +// if several polygons tile the region around a vertex, then exactly one of +// those polygons contains that vertex. Points and polylines still do not +// contain any vertices. +// +// - In the CLOSED model, all shapes contain their vertices (including points +// and polylines). +// +// Note that points other than vertices are never contained by polylines. +// If you want need this behavior, use S2ClosestEdgeQuery::IsDistanceLess() +// with a suitable distance threshold instead. +enum class S2VertexModel { OPEN, SEMI_OPEN, CLOSED }; + +// This class defines the options supported by S2ContainsPointQuery. +class S2ContainsPointQueryOptions { + public: + S2ContainsPointQueryOptions() {} + + // Convenience constructor that sets the vertex_model() option. + explicit S2ContainsPointQueryOptions(S2VertexModel vertex_model); + + // Controls whether shapes are considered to contain their vertices (see + // definitions above). By default the SEMI_OPEN model is used. + // + // DEFAULT: S2VertexModel::SEMI_OPEN + S2VertexModel vertex_model() const; + void set_vertex_model(S2VertexModel model); +}; +``` + +Here are the other methods defined by `S2ContainsPointQuery`: + +```c++ +// S2ContainsPointQuery determines whether one or more shapes in an +// S2ShapeIndex contain a given S2Point. The S2ShapeIndex may contain any +// number of points, polylines, and/or polygons (possibly overlapping). +// Shape boundaries may be modeled as OPEN, SEMI_OPEN, or CLOSED (this affects +// whether or not shapes are considered to contain their vertices). +// +// Example usage: +// S2ContainsPointQueryOptions options(S2VertexModel::CLOSED); +// return MakeS2ContainsPointQuery(&index, options).Contains(point); +// +// This class is not thread-safe. To use it in parallel, each thread should +// construct its own instance (this is not expensive). +// +// However, note that if you need to do a large number of point containment +// tests, it is more efficient to re-use the S2ContainsPointQuery object +// rather than constructing a new one each time. +template +class S2ContainsPointQuery { + public: + // Default constructor; requires Init() to be called. + S2ContainsPointQuery(); + + // Rather than calling this constructor, which requires specifying the + // IndexType template argument explicitly, the preferred idiom is to call + // MakeS2ContainsPointQuery() instead. For example: + // + // return MakeS2ContainsPointQuery(&index).Contains(p); + using Options = S2ContainsPointQueryOptions; + explicit S2ContainsPointQuery(const IndexType* index, + const Options& options = Options()); + + const IndexType& index() const; + const Options& options() const; + + // Equivalent to the two-argument constructor above. + void Init(const IndexType* index, const Options& options = Options()); + + // Returns true if any shape in the given index() contains the point "p" + // under the vertex model specified (OPEN, SEMI_OPEN, or CLOSED). + bool Contains(const S2Point& p); + + // Returns true if the given shape contains the point "p" under the vertex + // model specified (OPEN, SEMI_OPEN, or CLOSED). + // + // REQUIRES: "shape" belongs to index(). + bool ShapeContains(const S2Shape& shape, const S2Point& p); + + // Visits all shapes in the given index() that contain the given point "p", + // terminating early if the given ShapeVisitor function returns false (in + // which case VisitContainingShapes returns false as well). Each shape is + // visited at most once. + // + // Note that the API allows non-const access to the visited shapes. + using ShapeVisitor = std::function; + bool VisitContainingShapes(const S2Point& p, const ShapeVisitor& visitor); + + // Convenience function that returns all the shapes that contain the given + // point "p". + std::vector GetContainingShapes(const S2Point& p); + + // Visits all edges in the given index() that are incident to the point "p" + // (i.e., "p" is one of the edge endpoints), terminating early if the given + // EdgeVisitor function returns false (in which case VisitIncidentEdges + // returns false as well). Each edge is visited at most once. + using EdgeVisitor = std::function; + bool VisitIncidentEdges(const S2Point& p, const EdgeVisitor& visitor); +}; +``` + +## S2ShapeIndex Interface + +Most clients will not need to use the `S2ShapeIndex` interface directly, +since most operations can be accomplished by simply passing an +`S2ShapeIndex` to one of the classes mentioned in the overview +(`S2ClosestEdgeQuery`, `S2BooleanOperation`, etc). However the +`S2ShapeIndex` interface is briefly described here for the benefit of +clients that wish to implement their own algorithms. + +You can think of an `S2ShapeIndex` as a quad-tree that starts with the six +top-level faces of the `S2Cell` hierarchy and adaptively splits nodes that +intersect too many edges. However rather than storing the entire +quad-tree, instead `S2ShapeIndex` stores only the non-empty leaf nodes of +the quad-tree. This yields an ordered map whose keys are a set of +non-overlapping `S2CellIds`. For each `S2CellId`, the map stores the +subset of shapes and edges that intersect that `S2CellId`. + +Here is a visualization of an `S2ShapeIndex` for a polygon bounded by a +fractal curve: + + + +Notice that the `S2ShapeIndex` has entries for cells that are entirely in +the interior of the polygon. The existence of such cells makes point +location very fast. + +For each `S2CellId`, the index stores shapes and edges which intersect +that cell. An `S2ClippedShape` represents the intersection of a given +`S2Shape` with an `S2CellId`: + +```c++ +// S2ClippedShape represents the part of a shape that intersects an S2Cell. +// It consists of the set of edge ids that intersect that cell, and a boolean +// indicating whether the center of the cell is inside the shape (for shapes +// that have an interior). +// +// Note that the edges themselves are not clipped; we always use the original +// edges for intersection tests so that the results will be the same as the +// original shape. +class S2ClippedShape { + public: + // The shape id of the clipped shape. + int shape_id() const; + + // Return true if the center of the S2CellId is inside the shape. Returns + // false for shapes that do not have an interior. + bool contains_center() const; + + // The number of edges that intersect the S2CellId. + int num_edges() const; + + // Return the edge id of the given edge in this clipped shape. Edges are + // sorted in increasing order of edge id. + // + // REQUIRES: 0 <= i < num_edges() + int edge(int i) const; + + // Return true if the clipped shape contains the given edge id. + bool ContainsEdge(int id) const; +}; +``` + +An `S2ShapeIndexCell` represents the set of `S2ClippedShapes` that +intersect an `S2CellId': + +```c++ +// S2ShapeIndexCell stores the index contents for a given S2CellId. +// It consists of a set of clipped shapes. +class S2ShapeIndexCell { + public: + // Return the number of clipped shapes in this cell. + int num_clipped() const; + + // Return the clipped shape at the given index. Shapes are kept sorted in + // increasing order of shape id. + // + // REQUIRES: 0 <= i < num_clipped() + const S2ClippedShape& clipped(int i) const; + + // Return a pointer to the clipped shape corresponding to the given shape, + // or nullptr if the shape does not intersect this cell. + const S2ClippedShape* find_clipped(const S2Shape* shape) const; + const S2ClippedShape* find_clipped(int shape_id) const; + + // Convenience method that returns the total number of edges in all clipped + // shapes. + int num_edges() const; +}; +``` + +Finally, the `S2ShapeIndex` interface itself provides the following: + +* The ability to access the underlying `S2Shape` objects in the index. + Shapes can be retrieved by their *shape id*, an integer value that is + assigned when the shape is first added to the index (starting at zero + and increasing consecutively). + +* A random-access `Iterator` type that allows iterating and seeking over + `S2CellIds` in the index. For example, it can quickly determine the + `S2CellId` (if any) containing a given point. + +```c++ +class S2ShapeIndex { + public: + virtual ~S2ShapeIndex(); + + // Returns the number of distinct shape ids in the index. This is the same + // as the number of shapes provided that no shapes have ever been removed. + // (Shape ids are never reused.) + virtual int num_shape_ids() const = 0; + + // Returns a pointer to the shape with the given id, or nullptr if the shape + // has been removed from the index. + virtual S2Shape* shape(int id) const = 0; + + // Minimizes memory usage by requesting that any data structures that can be + // rebuilt should be discarded. This method invalidates all iterators. + // + // Like all non-const methods, this method is not thread-safe. + virtual void Minimize() = 0; + + // The possible relationships between a "target" cell and the cells of the + // S2ShapeIndex. If the target is an index cell or is contained by an index + // cell, it is "INDEXED". If the target is subdivided into one or more + // index cells, it is "SUBDIVIDED". Otherwise it is "DISJOINT". + enum CellRelation { + INDEXED, // Target is contained by an index cell + SUBDIVIDED, // Target is subdivided into one or more index cells + DISJOINT // Target does not intersect any index cells + }; + + // When passed to an Iterator constructor, specifies whether the iterator + // should be positioned at the beginning of the index (BEGIN), the end of + // the index (END), or arbitrarily (UNPOSITIONED). By default iterators are + // unpositioned, since this avoids an extra seek in this situation where one + // of the seek methods (such as Locate) is immediately called. + enum InitialPosition { BEGIN, END, UNPOSITIONED }; + + // A random access iterator that provides low-level access to the cells of + // the index. Cells are sorted in increasing order of S2CellId. + class Iterator { + public: + // Default constructor; must be followed by a call to Init(). + Iterator(); + + // Constructs an iterator positioned as specified. By default iterators + // are unpositioned, since this avoids an extra seek in this situation + // where one of the seek methods (such as Locate) is immediately called. + // + // If you want to position the iterator at the beginning, e.g. in order to + // loop through the entire index, do this instead: + // + // for (S2ShapeIndex::Iterator it(&index, S2ShapeIndex::BEGIN); + // !it.done(); it.Next()) { ... } + explicit Iterator(const S2ShapeIndex* index, + InitialPosition pos = UNPOSITIONED); + + // Initializes an iterator for the given S2ShapeIndex. This method may + // also be called in order to restore an iterator to a valid state after + // the underlying index has been updated (although it is usually easier + // just to declare a new iterator whenever required, since iterator + // construction is cheap). + void Init(const S2ShapeIndex* index, InitialPosition pos = UNPOSITIONED); + + // Iterators are copyable and moveable. + Iterator(const Iterator&); + Iterator& operator=(const Iterator&); + Iterator(Iterator&&); + Iterator& operator=(Iterator&&); + + // Returns the S2CellId of the current index cell. If done() is true, + // returns a value larger than any valid S2CellId (S2CellId::Sentinel()). + S2CellId id() const; + + // Returns the center point of the cell. + // REQUIRES: !done() + S2Point center() const; + + // Returns a reference to the contents of the current index cell. + // REQUIRES: !done() + const S2ShapeIndexCell& cell() const; + + // Returns true if the iterator is positioned past the last index cell. + bool done() const; + + // Positions the iterator at the first index cell (if any). + void Begin(); + + // Positions the iterator past the last index cell. + void Finish(); + + // Positions the iterator at the next index cell. + // REQUIRES: !done() + void Next(); + + // If the iterator is already positioned at the beginning, returns false. + // Otherwise positions the iterator at the previous entry and returns true. + bool Prev(); + + // Positions the iterator at the first cell with id() >= target, or at the + // end of the index if no such cell exists. + void Seek(S2CellId target); + + // Positions the iterator at the cell containing "target". If no such cell + // exists, returns false and leaves the iterator positioned arbitrarily. + // The returned index cell is guaranteed to contain all edges that might + // intersect the line segment between "target" and the cell center. + bool Locate(const S2Point& target); + + // Let T be the target S2CellId. If T is contained by some index cell I + // (including equality), this method positions the iterator at I and + // returns INDEXED. Otherwise if T contains one or more (smaller) index + // cells, it positions the iterator at the first such cell I and returns + // SUBDIVIDED. Otherwise it returns DISJOINT and leaves the iterator + // positioned arbitrarily. + CellRelation Locate(S2CellId target); + }; +``` + +Note that in addition to implementing a shared set of virtual methods, +`S2ShapeIndex` and its subtypes define an Iterator type with identical +methods. This makes it possible to support multiple `S2ShapeIndex` +subtypes either by accepting an instance of the base class: + + void DoSomething(const S2ShapeIndex& index) { + S2ShapeIndex::Iterator it(&index); + if (it.Locate(...)) { ... } + } + +or by accepting the `S2ShapeIndex` subtype as a template argument: + + template + void DoSomething(const IndexType& index) { + typename IndexType::Iterator it(&index); + if (it.Locate(...)) { ... } + } diff --git a/docs/favicon.ico b/docs/favicon.ico new file mode 100644 index 00000000..51699d93 Binary files /dev/null and b/docs/favicon.ico differ diff --git a/docs/favicons/android-chrome-192x192.png b/docs/favicons/android-chrome-192x192.png new file mode 100755 index 00000000..61c77570 Binary files /dev/null and b/docs/favicons/android-chrome-192x192.png differ diff --git a/docs/favicons/apple-touch-icon.png b/docs/favicons/apple-touch-icon.png new file mode 100755 index 00000000..3a4bed4d Binary files /dev/null and b/docs/favicons/apple-touch-icon.png differ diff --git a/docs/favicons/browserconfig.xml b/docs/favicons/browserconfig.xml new file mode 100755 index 00000000..3c3d111b --- /dev/null +++ b/docs/favicons/browserconfig.xml @@ -0,0 +1,9 @@ + + + + + + #ffffff + + + diff --git a/docs/favicons/favicon-16x16.png b/docs/favicons/favicon-16x16.png new file mode 100755 index 00000000..fc0654c0 Binary files /dev/null and b/docs/favicons/favicon-16x16.png differ diff --git a/docs/favicons/favicon-32x32.png b/docs/favicons/favicon-32x32.png new file mode 100755 index 00000000..76d9456f Binary files /dev/null and b/docs/favicons/favicon-32x32.png differ diff --git a/docs/favicons/favicon.ico b/docs/favicons/favicon.ico new file mode 100644 index 00000000..51699d93 Binary files /dev/null and b/docs/favicons/favicon.ico differ diff --git a/docs/favicons/manifest.json b/docs/favicons/manifest.json new file mode 100755 index 00000000..997639f6 --- /dev/null +++ b/docs/favicons/manifest.json @@ -0,0 +1,12 @@ +{ + "name": "grpc", + "icons": [ + { + "src": "\/android-chrome-192x192.png", + "sizes": "192x192", + "type": "image\/png" + } + ], + "theme_color": "#ffffff", + "display": "standalone" +} diff --git a/docs/favicons/mstile-150x150.png b/docs/favicons/mstile-150x150.png new file mode 100755 index 00000000..4ddc2efd Binary files /dev/null and b/docs/favicons/mstile-150x150.png differ diff --git a/docs/favicons/safari-pinned-tab.svg b/docs/favicons/safari-pinned-tab.svg new file mode 100755 index 00000000..b9af8d7c --- /dev/null +++ b/docs/favicons/safari-pinned-tab.svg @@ -0,0 +1,24 @@ + + + + +Created by potrace 1.11, written by Peter Selinger 2001-2013 + + + + + diff --git a/docs/getS2.md b/docs/getS2.md new file mode 100644 index 00000000..5065e949 --- /dev/null +++ b/docs/getS2.md @@ -0,0 +1,22 @@ +--- +title: Get S2 +--- + +Code for S2 is located in the following repositories: + +Reference Implementations: + +* C++ S2 Library + +Ported Subsets: + +* Java + S2 Library +* Go S2 Library +* Python + S2 Library + +Clone any of these repositories to get started working with S2! You may wish +to fork these repositories if you wish to make contributions (where applicable). + +Consult the [Platforms Guide](about/platforms) for installation instructions. diff --git a/docs/img/corner_10.gif b/docs/img/corner_10.gif new file mode 100644 index 00000000..ea19c14e Binary files /dev/null and b/docs/img/corner_10.gif differ diff --git a/docs/img/corner_1000.gif b/docs/img/corner_1000.gif new file mode 100644 index 00000000..23d4c64d Binary files /dev/null and b/docs/img/corner_1000.gif differ diff --git a/docs/img/corner_20.gif b/docs/img/corner_20.gif new file mode 100644 index 00000000..f7127cd3 Binary files /dev/null and b/docs/img/corner_20.gif differ diff --git a/docs/img/corner_200.gif b/docs/img/corner_200.gif new file mode 100644 index 00000000..f94a39da Binary files /dev/null and b/docs/img/corner_200.gif differ diff --git a/docs/img/corner_50.gif b/docs/img/corner_50.gif new file mode 100644 index 00000000..ca6cb4ab Binary files /dev/null and b/docs/img/corner_50.gif differ diff --git a/docs/img/edge_4.gif b/docs/img/edge_4.gif new file mode 100644 index 00000000..a6bbed50 Binary files /dev/null and b/docs/img/edge_4.gif differ diff --git a/docs/img/edge_50.gif b/docs/img/edge_50.gif new file mode 100644 index 00000000..965a3a6f Binary files /dev/null and b/docs/img/edge_50.gif differ diff --git a/docs/img/edge_500.gif b/docs/img/edge_500.gif new file mode 100644 index 00000000..1707798f Binary files /dev/null and b/docs/img/edge_500.gif differ diff --git a/docs/img/edge_6.gif b/docs/img/edge_6.gif new file mode 100644 index 00000000..60696afa Binary files /dev/null and b/docs/img/edge_6.gif differ diff --git a/docs/img/face0.jpg b/docs/img/face0.jpg new file mode 100644 index 00000000..fd929e33 Binary files /dev/null and b/docs/img/face0.jpg differ diff --git a/docs/img/face0_disp.jpg b/docs/img/face0_disp.jpg new file mode 100644 index 00000000..fd3bf8e5 Binary files /dev/null and b/docs/img/face0_disp.jpg differ diff --git a/docs/img/face1.jpg b/docs/img/face1.jpg new file mode 100644 index 00000000..ff8e5a03 Binary files /dev/null and b/docs/img/face1.jpg differ diff --git a/docs/img/face1_disp.jpg b/docs/img/face1_disp.jpg new file mode 100644 index 00000000..f2d26b24 Binary files /dev/null and b/docs/img/face1_disp.jpg differ diff --git a/docs/img/face2.jpg b/docs/img/face2.jpg new file mode 100644 index 00000000..6ab8272a Binary files /dev/null and b/docs/img/face2.jpg differ diff --git a/docs/img/face2_disp.jpg b/docs/img/face2_disp.jpg new file mode 100644 index 00000000..306cff46 Binary files /dev/null and b/docs/img/face2_disp.jpg differ diff --git a/docs/img/face3.jpg b/docs/img/face3.jpg new file mode 100644 index 00000000..474bd558 Binary files /dev/null and b/docs/img/face3.jpg differ diff --git a/docs/img/face3_disp.jpg b/docs/img/face3_disp.jpg new file mode 100644 index 00000000..3ff2164e Binary files /dev/null and b/docs/img/face3_disp.jpg differ diff --git a/docs/img/face4.jpg b/docs/img/face4.jpg new file mode 100644 index 00000000..884e2fe9 Binary files /dev/null and b/docs/img/face4.jpg differ diff --git a/docs/img/face4_disp.jpg b/docs/img/face4_disp.jpg new file mode 100644 index 00000000..56b6da0c Binary files /dev/null and b/docs/img/face4_disp.jpg differ diff --git a/docs/img/face5.jpg b/docs/img/face5.jpg new file mode 100644 index 00000000..d85fb47e Binary files /dev/null and b/docs/img/face5.jpg differ diff --git a/docs/img/face5_disp.jpg b/docs/img/face5_disp.jpg new file mode 100644 index 00000000..54411c21 Binary files /dev/null and b/docs/img/face5_disp.jpg differ diff --git a/docs/img/florida1.gif b/docs/img/florida1.gif new file mode 100644 index 00000000..9885e8ce Binary files /dev/null and b/docs/img/florida1.gif differ diff --git a/docs/img/florida2.gif b/docs/img/florida2.gif new file mode 100644 index 00000000..d0820062 Binary files /dev/null and b/docs/img/florida2.gif differ diff --git a/docs/img/hawaii.gif b/docs/img/hawaii.gif new file mode 100644 index 00000000..298182d5 Binary files /dev/null and b/docs/img/hawaii.gif differ diff --git a/docs/img/hilbert-cell-center.gif b/docs/img/hilbert-cell-center.gif new file mode 100644 index 00000000..f41007f3 Binary files /dev/null and b/docs/img/hilbert-cell-center.gif differ diff --git a/docs/img/hilbert-figure.gif b/docs/img/hilbert-figure.gif new file mode 100644 index 00000000..2fb83ad4 Binary files /dev/null and b/docs/img/hilbert-figure.gif differ diff --git a/docs/img/kirkland_1.gif b/docs/img/kirkland_1.gif new file mode 100644 index 00000000..24e1811f Binary files /dev/null and b/docs/img/kirkland_1.gif differ diff --git a/docs/img/kirkland_10.gif b/docs/img/kirkland_10.gif new file mode 100644 index 00000000..aaf10aa6 Binary files /dev/null and b/docs/img/kirkland_10.gif differ diff --git a/docs/img/kirkland_2.gif b/docs/img/kirkland_2.gif new file mode 100644 index 00000000..629a31e5 Binary files /dev/null and b/docs/img/kirkland_2.gif differ diff --git a/docs/img/kirkland_20.gif b/docs/img/kirkland_20.gif new file mode 100644 index 00000000..b0cb6663 Binary files /dev/null and b/docs/img/kirkland_20.gif differ diff --git a/docs/img/kirkland_3.gif b/docs/img/kirkland_3.gif new file mode 100644 index 00000000..fa810824 Binary files /dev/null and b/docs/img/kirkland_3.gif differ diff --git a/docs/img/kirkland_4.gif b/docs/img/kirkland_4.gif new file mode 100644 index 00000000..ff3d0dcd Binary files /dev/null and b/docs/img/kirkland_4.gif differ diff --git a/docs/img/kirkland_5.gif b/docs/img/kirkland_5.gif new file mode 100644 index 00000000..e3e1ab3b Binary files /dev/null and b/docs/img/kirkland_5.gif differ diff --git a/docs/img/kirkland_50.gif b/docs/img/kirkland_50.gif new file mode 100644 index 00000000..3ecb889b Binary files /dev/null and b/docs/img/kirkland_50.gif differ diff --git a/docs/img/kirkland_500.gif b/docs/img/kirkland_500.gif new file mode 100644 index 00000000..11ca930a Binary files /dev/null and b/docs/img/kirkland_500.gif differ diff --git a/docs/img/kirkland_6.gif b/docs/img/kirkland_6.gif new file mode 100644 index 00000000..e5a036c5 Binary files /dev/null and b/docs/img/kirkland_6.gif differ diff --git a/docs/img/polar_100.gif b/docs/img/polar_100.gif new file mode 100644 index 00000000..eba23b05 Binary files /dev/null and b/docs/img/polar_100.gif differ diff --git a/docs/img/polar_20.gif b/docs/img/polar_20.gif new file mode 100644 index 00000000..41411c1b Binary files /dev/null and b/docs/img/polar_20.gif differ diff --git a/docs/img/polar_500.gif b/docs/img/polar_500.gif new file mode 100644 index 00000000..5431bccc Binary files /dev/null and b/docs/img/polar_500.gif differ diff --git a/docs/img/polar_8.gif b/docs/img/polar_8.gif new file mode 100644 index 00000000..35d6881d Binary files /dev/null and b/docs/img/polar_8.gif differ diff --git a/docs/img/s2cell_global.jpg b/docs/img/s2cell_global.jpg new file mode 100644 index 00000000..ff9e3e2f Binary files /dev/null and b/docs/img/s2cell_global.jpg differ diff --git a/docs/img/s2curve-large.gif b/docs/img/s2curve-large.gif new file mode 100644 index 00000000..937871c1 Binary files /dev/null and b/docs/img/s2curve-large.gif differ diff --git a/docs/img/s2curve-medium.gif b/docs/img/s2curve-medium.gif new file mode 100644 index 00000000..f7a2d569 Binary files /dev/null and b/docs/img/s2curve-medium.gif differ diff --git a/docs/img/s2curve-small.gif b/docs/img/s2curve-small.gif new file mode 100644 index 00000000..f75ca234 Binary files /dev/null and b/docs/img/s2curve-small.gif differ diff --git a/docs/img/s2curve-tiny.gif b/docs/img/s2curve-tiny.gif new file mode 100644 index 00000000..0c1fd7c8 Binary files /dev/null and b/docs/img/s2curve-tiny.gif differ diff --git a/docs/img/s2hierarchy.gif b/docs/img/s2hierarchy.gif new file mode 100644 index 00000000..711aa35a Binary files /dev/null and b/docs/img/s2hierarchy.gif differ diff --git a/docs/img/s2shapeindex_example.png b/docs/img/s2shapeindex_example.png new file mode 100644 index 00000000..b7cc8e31 Binary files /dev/null and b/docs/img/s2shapeindex_example.png differ diff --git a/docs/img/washington_10.gif b/docs/img/washington_10.gif new file mode 100644 index 00000000..77794e9e Binary files /dev/null and b/docs/img/washington_10.gif differ diff --git a/docs/img/washington_100.gif b/docs/img/washington_100.gif new file mode 100644 index 00000000..bcf0e19f Binary files /dev/null and b/docs/img/washington_100.gif differ diff --git a/docs/img/washington_20.gif b/docs/img/washington_20.gif new file mode 100644 index 00000000..22aedaaa Binary files /dev/null and b/docs/img/washington_20.gif differ diff --git a/docs/img/xyz_to_uv_to_st.gif b/docs/img/xyz_to_uv_to_st.gif new file mode 100644 index 00000000..7ec9fa45 Binary files /dev/null and b/docs/img/xyz_to_uv_to_st.gif differ diff --git a/docs/index.md b/docs/index.md new file mode 100644 index 00000000..421f8654 --- /dev/null +++ b/docs/index.md @@ -0,0 +1,72 @@ +--- +title: S2 Geometry +--- + + +
Welcome to the S2 Geometry Library!
+
+ +![](devguide/img/s2hierarchy.gif) + +A unique feature of the S2 library is that unlike traditional geographic +information systems, which represent data as flat two-dimensional projections +(similar to an atlas), the S2 library represents all data on a three-dimensional +sphere (similar to a globe). This makes it possible to build a worldwide +geographic database with no seams or singularities, using a single coordinate +system, and with low distortion everywhere compared to the true shape of the +Earth. While the Earth is not quite spherical, it is much closer to being a +sphere than it is to being flat! + +If you want to learn more about the library, start by reading the +[Overview](about/overview) and +[QuickStart document](devguide/cpp/quickstart), then read the +introduction to the [basic types](devguide/basic_types) + +** +Get S2 on GitHub** + +## S2 Features + +Notable features of the library include: + +* Flexible support for spatial indexing, including the ability + to approximate arbitrary regions as collections of discrete + S2 cells. This feature makes it easy to build large distributed + spatial indexes. +* Fast in-memory spatial indexing of collections of points, + polylines, and polygons. +* Robust constructive operations (such as intersection, union, and + simplification) and boolean predicates (such as testing for + containment). +* Efficient query operations for finding nearby objects, measuring + distances, computing centroids, etc. +* A flexible and robust implementation of snap rounding (a geometric + technique that allows operations to be implemented 100% robustly + while using small and fast coordinate representations). +* A collection of efficient yet exact mathematical predicates for + testing relationships among geometric primitives. +* Extensive testing on Google's vast collection of geographic data. +* Flexible Apache 2.0 license. + +## Table of Contents + +* About S2 + + * [Overview](about/overview) + * [Platforms Guide](about/platforms) + +* Tutorials + + * [Quick Start](devguide/cpp/quickstart) + +* Developer Guides + + * [Basic Types](devguide/basic_types) + * [Finding Nearby Edges](devguide/s2closestedgequery) + * [S2Cell Hierarchy](devguide/s2cell_hierarchy) + * [Covering Examples](devguide/examples/coverings) + +* Resources + + * The [S2 Earthcube](resources/earthcube) + * [S2Cell Statistics](resources/s2cell_statistics) diff --git a/docs/resources/earthcube.md b/docs/resources/earthcube.md new file mode 100644 index 00000000..b1110366 --- /dev/null +++ b/docs/resources/earthcube.md @@ -0,0 +1,51 @@ +--- +title: Earth Cube +--- + +The S2Cell hierarchy projects the Earth onto six top-level "face cells", +which may be unwrapped and flattened into a stairstep arrangement. Below +is a map of the six faces: + +![](/devguide/img/s2cell_global.jpg) + +A single space-filling curve covers the entire globe, snaking from one +face to the next. Note the traversal order in the odd-numbered faces is +the mirror of the order in even numbered faces. + +~~~~ + ^----> <----^ + | 4 | 5 | + | | | + x V -----> + + ^----> <----^ + | 2 | 3 | + | | | + x V -----> + + ^----> <----^ + | 0 | 1 | + | | | + x v -----> +~~~~ + +You can print out the following six images on a good color printer, glue +them onto sturdy cardboard squares and tape them together to form a +cube, to make a model S2Cell-Hierarchy Globe, which will let you +visually locate any cell at level 5 or higher. + +CAVEAT: While these images are conceptually useful, the satellite imagery has +not been transformed correctly, so don't rely on these images to determine +which geographic features are contained by a given `S2CellId`. + +[![](img/face0_disp.jpg)](img/face0.jpg) + +[![](img/face1_disp.jpg)](img/face1.jpg) + +[![](img/face2_disp.jpg)](img/face2.jpg) + +[![](img/face3_disp.jpg)](img/face3.jpg) + +[![](img/face4_disp.jpg)](img/face4.jpg) + +[![](img/face5_disp.jpg)](img/face5.jpg) diff --git a/docs/resources/img/face0.jpg b/docs/resources/img/face0.jpg new file mode 100644 index 00000000..fd929e33 Binary files /dev/null and b/docs/resources/img/face0.jpg differ diff --git a/docs/resources/img/face0_disp.jpg b/docs/resources/img/face0_disp.jpg new file mode 100644 index 00000000..fd3bf8e5 Binary files /dev/null and b/docs/resources/img/face0_disp.jpg differ diff --git a/docs/resources/img/face1.jpg b/docs/resources/img/face1.jpg new file mode 100644 index 00000000..ff8e5a03 Binary files /dev/null and b/docs/resources/img/face1.jpg differ diff --git a/docs/resources/img/face1_disp.jpg b/docs/resources/img/face1_disp.jpg new file mode 100644 index 00000000..f2d26b24 Binary files /dev/null and b/docs/resources/img/face1_disp.jpg differ diff --git a/docs/resources/img/face2.jpg b/docs/resources/img/face2.jpg new file mode 100644 index 00000000..6ab8272a Binary files /dev/null and b/docs/resources/img/face2.jpg differ diff --git a/docs/resources/img/face2_disp.jpg b/docs/resources/img/face2_disp.jpg new file mode 100644 index 00000000..306cff46 Binary files /dev/null and b/docs/resources/img/face2_disp.jpg differ diff --git a/docs/resources/img/face3.jpg b/docs/resources/img/face3.jpg new file mode 100644 index 00000000..474bd558 Binary files /dev/null and b/docs/resources/img/face3.jpg differ diff --git a/docs/resources/img/face3_disp.jpg b/docs/resources/img/face3_disp.jpg new file mode 100644 index 00000000..3ff2164e Binary files /dev/null and b/docs/resources/img/face3_disp.jpg differ diff --git a/docs/resources/img/face4.jpg b/docs/resources/img/face4.jpg new file mode 100644 index 00000000..884e2fe9 Binary files /dev/null and b/docs/resources/img/face4.jpg differ diff --git a/docs/resources/img/face4_disp.jpg b/docs/resources/img/face4_disp.jpg new file mode 100644 index 00000000..56b6da0c Binary files /dev/null and b/docs/resources/img/face4_disp.jpg differ diff --git a/docs/resources/img/face5.jpg b/docs/resources/img/face5.jpg new file mode 100644 index 00000000..d85fb47e Binary files /dev/null and b/docs/resources/img/face5.jpg differ diff --git a/docs/resources/img/face5_disp.jpg b/docs/resources/img/face5_disp.jpg new file mode 100644 index 00000000..54411c21 Binary files /dev/null and b/docs/resources/img/face5_disp.jpg differ diff --git a/docs/resources/img/s2cell_global.jpg b/docs/resources/img/s2cell_global.jpg new file mode 100644 index 00000000..ff9e3e2f Binary files /dev/null and b/docs/resources/img/s2cell_global.jpg differ diff --git a/docs/resources/index.md b/docs/resources/index.md new file mode 100644 index 00000000..abbc2539 --- /dev/null +++ b/docs/resources/index.md @@ -0,0 +1,6 @@ +--- +title: Resources +--- + +* [S2 Cell Statistics](s2cell_statistics) +* [S2 Earthcube](earthcube) diff --git a/docs/resources/s2cell_statistics.md b/docs/resources/s2cell_statistics.md new file mode 100644 index 00000000..a294ae56 --- /dev/null +++ b/docs/resources/s2cell_statistics.md @@ -0,0 +1,44 @@ +--- +title: S2 Cell Statistics +--- + +These are the approximate ranges of areas of +[S2 cells](/about/overview.md#s2cell-hierarchy) +at each level. The average size is that returned by `S2Cell.AverageArea()`, +guaranteed to be within a factor of 1.5 of the high and low end of the range. +Additionally, two random locations were chosen to provide example cell edge +lengths. + +| **level** | **min area** | **max area** | **average area** | **units** | | **Random cell 1 (UK) min edge length** | **Random cell 1 (UK) max edge length** | | **Random cell 2 (US) min edge length** | **Random cell 2 (US) max edge length** | **Number of cells** | +| --- | ---: | ---: | ---: | ---: | --- | ---: | ---: | --- | ---: | ---: | ---: | +| 00 | 85011012.19 | 85011012.19 | 85011012.19 | km2 | | 7842 km | 7842 km | | 7842 km | 7842 km | 6 | +| 01 | 21252753.05 | 21252753.05 | 21252753.05 | km2 | | 3921 km | 5004 km | | 3921 km | 5004 km | 24 | +| 02 | 4919708.23 | 6026521.16 | 5313188.26 | km2 | | 1825 km | 2489 km | | 1825 km | 2489 km | 96 | +| 03 | 1055377.48 | 1646455.50 | 1328297.07 | km2 | | 840 km | 1167 km | | 1130 km | 1310 km | 384 | +| 04 | 231564.06 | 413918.15 | 332074.27 | km2 | | 432 km | 609 km | | 579 km | 636 km | 1536 | +| 05 | 53798.67 | 104297.91 | 83018.57 | km2 | | 210 km | 298 km | | 287 km | 315 km | 6K | +| 06 | 12948.81 | 26113.30 | 20754.64 | km2 | | 108 km | 151 km | | 143 km | 156 km | 24K | +| 07 | 3175.44 | 6529.09 | 5188.66 | km2 | | 54 km | 76 km | | 72 km | 78 km | 98K | +| 08 | 786.20 | 1632.45 | 1297.17 | km2 | | 27 km | 38 km | | 36 km | 39 km | 393K | +| 09 | 195.59 | 408.12 | 324.29 | km2 | | 14 km | 19 km | | 18 km | 20 km | 1573K | +| 10 | 48.78 | 102.03 | 81.07 | km2 | | 7 km | 9 km | | 9 km | 10 km | 6M | +| 11 | 12.18 | 25.51 | 20.27 | km2 | | 3 km | 5 km | | 4 km | 5 km | 25M | +| 12 | 3.04 | 6.38 | 5.07 | km2 | | 1699 m | 2 km | | 2 km | 2 km | 100M | +| 13 | 0.76 | 1.59 | 1.27 | km2 | | 850 m | 1185 m | | 1123 m | 1225 m | 402M | +| 14 | 0.19 | 0.40 | 0.32 | km2 | | 425 m | 593 m | | 562 m | 613 m | 1610M | +| 15 | 47520.30 | 99638.93 | 79172.67 | m2 | | 212 m | 296 m | | 281 m | 306 m | 6B | +| 16 | 11880.08 | 24909.73 | 19793.17 | m2 | | 106 m | 148 m | | 140 m | 153 m | 25B | +| 17 | 2970.02 | 6227.43 | 4948.29 | m2 | | 53 m | 74 m | | 70 m | 77 m | 103B | +| 18 | 742.50 | 1556.86 | 1237.07 | m2 | | 27 m | 37 m | | 35 m | 38 m | 412B | +| 19 | 185.63 | 389.21 | 309.27 | m2 | | 13 m | 19 m | | 18 m | 19 m | 1649B | +| 20 | 46.41 | 97.30 | 77.32 | m2 | | 7 m | 9 m | | 9 m | 10 m | 7T | +| 21 | 11.60 | 24.33 | 19.33 | m2 | | 3 m | 5 m | | 4 m | 5 m | 26T | +| 22 | 2.90 | 6.08 | 4.83 | m2 | | 166 cm | 2 m | | 2 m | 2 m | 105T | +| 23 | 0.73 | 1.52 | 1.21 | m2 | | 83 cm | 116 cm | | 110 cm | 120 cm | 422T | +| 24 | 0.18 | 0.38 | 0.30 | m2 | | 41 cm | 58 cm | | 55 cm | 60 cm | 1689T | +| 25 | 453.19 | 950.23 | 755.05 | cm2 | | 21 cm | 29 cm | | 27 cm | 30 cm | 7e15 | +| 26 | 113.30 | 237.56 | 188.76 | cm2 | | 10 cm | 14 cm | | 14 cm | 15 cm | 27e15 | +| 27 | 28.32 | 59.39 | 47.19 | cm2 | | 5 cm | 7 cm | | 7 cm | 7 cm | 108e15 | +| 28 | 7.08 | 14.85 | 11.80 | cm2 | | 2 cm | 4 cm | | 3 cm | 4 cm | 432e15 | +| 29 | 1.77 | 3.71 | 2.95 | cm2 | | 12 mm | 18 mm | | 17 mm | 18 mm | 1729e15 | +| 30 | 0.44 | 0.93 | 0.74 | cm2 | | 6 mm | 9 mm | | 8 mm | 9 mm | 7e18 |