diff options
Diffstat (limited to 'vendor/github.com/golang/geo/s2/regioncoverer.go')
| -rw-r--r-- | vendor/github.com/golang/geo/s2/regioncoverer.go | 477 |
1 files changed, 477 insertions, 0 deletions
diff --git a/vendor/github.com/golang/geo/s2/regioncoverer.go b/vendor/github.com/golang/geo/s2/regioncoverer.go new file mode 100644 index 000000000..476e58559 --- /dev/null +++ b/vendor/github.com/golang/geo/s2/regioncoverer.go @@ -0,0 +1,477 @@ +// Copyright 2015 Google Inc. All rights reserved. +// +// 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. + +package s2 + +import ( + "container/heap" +) + +// RegionCoverer allows arbitrary regions to be approximated as unions of cells (CellUnion). +// This is useful for implementing various sorts of search and precomputation operations. +// +// Typical usage: +// +// rc := &s2.RegionCoverer{MaxLevel: 30, MaxCells: 5} +// r := s2.Region(CapFromCenterArea(center, area)) +// covering := rc.Covering(r) +// +// This yields a CellUnion of at most 5 cells that is guaranteed to cover the +// given region (a disc-shaped region on the sphere). +// +// For covering, only cells where (level - MinLevel) is a multiple of LevelMod will be used. +// This effectively allows the branching factor of the S2 CellID 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. +// +// Note the following: +// +// - MinLevel takes priority over MaxCells, i.e. cells below the given level will +// never be used even if this causes a large number of cells to be returned. +// +// - For any setting of MaxCells, up to 6 cells may be returned if that +// is the minimum number of cells required (e.g. if the region intersects +// all six face cells). Up to 3 cells may be returned even for very tiny +// convex regions if they happen to be located at the intersection of +// three cube faces. +// +// - For any setting of MaxCells, an arbitrary number of cells may be +// returned if MinLevel is too high for the region being approximated. +// +// - If MaxCells is less than 4, the area of the covering may be +// arbitrarily large compared to the area of the original region even if +// the region is convex (e.g. a Cap or Rect). +// +// 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 MaxCells is a limit on how much work is done exploring the +// possible covering as well as a limit on the final output size. +// +// Because it is an approximation algorithm, one should not rely on the +// stability of the output. In particular, the output of the covering algorithm +// may change across different versions of the library. +// +// One can also generate interior coverings, which are sets of cells which +// are entirely contained within a region. Interior coverings can be +// empty, even for non-empty regions, if there are no cells that satisfy +// the provided constraints and are contained by the region. Note that for +// performance reasons, it is wise to specify a MaxLevel when computing +// interior coverings - otherwise for regions with small or zero area, the +// algorithm may spend a lot of time subdividing cells all the way to leaf +// level to try to find contained cells. +type RegionCoverer struct { + MinLevel int // the minimum cell level to be used. + MaxLevel int // the maximum cell level to be used. + LevelMod int // the LevelMod to be used. + MaxCells int // the maximum desired number of cells in the approximation. +} + +type coverer struct { + minLevel int // the minimum cell level to be used. + maxLevel int // the maximum cell level to be used. + levelMod int // the LevelMod to be used. + maxCells int // the maximum desired number of cells in the approximation. + region Region + result CellUnion + pq priorityQueue + interiorCovering bool +} + +type candidate struct { + cell Cell + terminal bool // Cell should not be expanded further. + numChildren int // Number of children that intersect the region. + children []*candidate // Actual size may be 0, 4, 16, or 64 elements. + priority int // Priority of the candidate. +} + +type priorityQueue []*candidate + +func (pq priorityQueue) Len() int { + return len(pq) +} + +func (pq priorityQueue) Less(i, j int) bool { + // We want Pop to give us the highest, not lowest, priority so we use greater than here. + return pq[i].priority > pq[j].priority +} + +func (pq priorityQueue) Swap(i, j int) { + pq[i], pq[j] = pq[j], pq[i] +} + +func (pq *priorityQueue) Push(x interface{}) { + item := x.(*candidate) + *pq = append(*pq, item) +} + +func (pq *priorityQueue) Pop() interface{} { + item := (*pq)[len(*pq)-1] + *pq = (*pq)[:len(*pq)-1] + return item +} + +func (pq *priorityQueue) Reset() { + *pq = (*pq)[:0] +} + +// newCandidate returns a new candidate with no children if the cell intersects the given region. +// The candidate is marked as terminal if it should not be expanded further. +func (c *coverer) newCandidate(cell Cell) *candidate { + if !c.region.IntersectsCell(cell) { + return nil + } + cand := &candidate{cell: cell} + level := int(cell.level) + if level >= c.minLevel { + if c.interiorCovering { + if c.region.ContainsCell(cell) { + cand.terminal = true + } else if level+c.levelMod > c.maxLevel { + return nil + } + } else if level+c.levelMod > c.maxLevel || c.region.ContainsCell(cell) { + cand.terminal = true + } + } + return cand +} + +// expandChildren populates the children of the candidate by expanding the given number of +// levels from the given cell. Returns the number of children that were marked "terminal". +func (c *coverer) expandChildren(cand *candidate, cell Cell, numLevels int) int { + numLevels-- + var numTerminals int + last := cell.id.ChildEnd() + for ci := cell.id.ChildBegin(); ci != last; ci = ci.Next() { + childCell := CellFromCellID(ci) + if numLevels > 0 { + if c.region.IntersectsCell(childCell) { + numTerminals += c.expandChildren(cand, childCell, numLevels) + } + continue + } + if child := c.newCandidate(childCell); child != nil { + cand.children = append(cand.children, child) + cand.numChildren++ + if child.terminal { + numTerminals++ + } + } + } + return numTerminals +} + +// addCandidate adds the given candidate to the result if it is marked as "terminal", +// otherwise expands its children and inserts it into the priority queue. +// Passing an argument of nil does nothing. +func (c *coverer) addCandidate(cand *candidate) { + if cand == nil { + return + } + + if cand.terminal { + c.result = append(c.result, cand.cell.id) + return + } + + // Expand one level at a time until we hit minLevel to ensure that we don't skip over it. + numLevels := c.levelMod + level := int(cand.cell.level) + if level < c.minLevel { + numLevels = 1 + } + + numTerminals := c.expandChildren(cand, cand.cell, numLevels) + maxChildrenShift := uint(2 * c.levelMod) + if cand.numChildren == 0 { + return + } else if !c.interiorCovering && numTerminals == 1<<maxChildrenShift && level >= c.minLevel { + // Optimization: add the parent cell rather than all of its children. + // We can't do this for interior coverings, since the children just + // intersect the region, but may not be contained by it - we need to + // subdivide them further. + cand.terminal = true + c.addCandidate(cand) + } else { + // We negate the priority so that smaller absolute priorities are returned + // first. The heuristic is designed to refine the largest cells first, + // since those are where we have the largest potential gain. Among cells + // of the same size, we prefer the cells with the fewest children. + // Finally, among cells with equal numbers of children we prefer those + // with the smallest number of children that cannot be refined further. + cand.priority = -(((level<<maxChildrenShift)+cand.numChildren)<<maxChildrenShift + numTerminals) + heap.Push(&c.pq, cand) + } +} + +// adjustLevel returns the reduced "level" so that it satisfies levelMod. Levels smaller than minLevel +// are not affected (since cells at these levels are eventually expanded). +func (c *coverer) adjustLevel(level int) int { + if c.levelMod > 1 && level > c.minLevel { + level -= (level - c.minLevel) % c.levelMod + } + return level +} + +// adjustCellLevels ensures that all cells with level > minLevel also satisfy levelMod, +// by replacing them with an ancestor if necessary. Cell levels smaller +// than minLevel are not modified (see AdjustLevel). The output is +// then normalized to ensure that no redundant cells are present. +func (c *coverer) adjustCellLevels(cells *CellUnion) { + if c.levelMod == 1 { + return + } + + var out int + for _, ci := range *cells { + level := ci.Level() + newLevel := c.adjustLevel(level) + if newLevel != level { + ci = ci.Parent(newLevel) + } + if out > 0 && (*cells)[out-1].Contains(ci) { + continue + } + for out > 0 && ci.Contains((*cells)[out-1]) { + out-- + } + (*cells)[out] = ci + out++ + } + *cells = (*cells)[:out] +} + +// initialCandidates computes a set of initial candidates that cover the given region. +func (c *coverer) initialCandidates() { + // Optimization: start with a small (usually 4 cell) covering of the region's bounding cap. + temp := &RegionCoverer{MaxLevel: c.maxLevel, LevelMod: 1, MaxCells: minInt(4, c.maxCells)} + + cells := temp.FastCovering(c.region) + c.adjustCellLevels(&cells) + for _, ci := range cells { + c.addCandidate(c.newCandidate(CellFromCellID(ci))) + } +} + +// coveringInternal generates a covering and stores it in result. +// Strategy: Start with the 6 faces of the cube. Discard any +// that do not intersect the shape. Then repeatedly choose the +// largest cell that intersects the shape and subdivide it. +// +// result contains the cells that will be part of the output, while pq +// contains cells that we may still subdivide further. Cells that are +// entirely contained within the region are immediately added to the output, +// while cells that do not intersect the region are immediately discarded. +// Therefore pq only contains cells that partially intersect the region. +// Candidates are prioritized first according to cell size (larger cells +// first), then by the number of intersecting children they have (fewest +// children first), and then by the number of fully contained children +// (fewest children first). +func (c *coverer) coveringInternal(region Region) { + c.region = region + + c.initialCandidates() + for c.pq.Len() > 0 && (!c.interiorCovering || len(c.result) < c.maxCells) { + cand := heap.Pop(&c.pq).(*candidate) + + // For interior covering we keep subdividing no matter how many children + // candidate has. If we reach MaxCells before expanding all children, + // we will just use some of them. + // For exterior covering we cannot do this, because result has to cover the + // whole region, so all children have to be used. + // candidate.numChildren == 1 case takes care of the situation when we + // already have more than MaxCells in result (minLevel is too high). + // Subdividing of the candidate with one child does no harm in this case. + if c.interiorCovering || int(cand.cell.level) < c.minLevel || cand.numChildren == 1 || len(c.result)+c.pq.Len()+cand.numChildren <= c.maxCells { + for _, child := range cand.children { + if !c.interiorCovering || len(c.result) < c.maxCells { + c.addCandidate(child) + } + } + } else { + cand.terminal = true + c.addCandidate(cand) + } + } + c.pq.Reset() + c.region = nil +} + +// newCoverer returns an instance of coverer. +func (rc *RegionCoverer) newCoverer() *coverer { + return &coverer{ + minLevel: maxInt(0, minInt(maxLevel, rc.MinLevel)), + maxLevel: maxInt(0, minInt(maxLevel, rc.MaxLevel)), + levelMod: maxInt(1, minInt(3, rc.LevelMod)), + maxCells: rc.MaxCells, + } +} + +// Covering returns a CellUnion that covers the given region and satisfies the various restrictions. +func (rc *RegionCoverer) Covering(region Region) CellUnion { + covering := rc.CellUnion(region) + covering.Denormalize(maxInt(0, minInt(maxLevel, rc.MinLevel)), maxInt(1, minInt(3, rc.LevelMod))) + return covering +} + +// InteriorCovering returns a CellUnion that is contained within the given region and satisfies the various restrictions. +func (rc *RegionCoverer) InteriorCovering(region Region) CellUnion { + intCovering := rc.InteriorCellUnion(region) + intCovering.Denormalize(maxInt(0, minInt(maxLevel, rc.MinLevel)), maxInt(1, minInt(3, rc.LevelMod))) + return intCovering +} + +// CellUnion returns a normalized CellUnion that covers the given region and +// satisfies the restrictions except for minLevel and levelMod. These criteria +// cannot be satisfied using a cell union because cell unions are +// automatically normalized by replacing four child cells with their parent +// whenever possible. (Note that the list of cell ids passed to the CellUnion +// constructor does in fact satisfy all the given restrictions.) +func (rc *RegionCoverer) CellUnion(region Region) CellUnion { + c := rc.newCoverer() + c.coveringInternal(region) + cu := c.result + cu.Normalize() + return cu +} + +// InteriorCellUnion returns a normalized CellUnion that is contained within the given region and +// satisfies the restrictions except for minLevel and levelMod. These criteria +// cannot be satisfied using a cell union because cell unions are +// automatically normalized by replacing four child cells with their parent +// whenever possible. (Note that the list of cell ids passed to the CellUnion +// constructor does in fact satisfy all the given restrictions.) +func (rc *RegionCoverer) InteriorCellUnion(region Region) CellUnion { + c := rc.newCoverer() + c.interiorCovering = true + c.coveringInternal(region) + cu := c.result + cu.Normalize() + return cu +} + +// FastCovering returns a CellUnion that covers the given region similar to Covering, +// except that this method is much faster and the coverings are not as tight. +// All of the usual parameters are respected (MaxCells, MinLevel, MaxLevel, and LevelMod), +// except that the implementation makes no attempt to take advantage of large values of +// MaxCells. (A small number of cells will always be returned.) +// +// This function is useful as a starting point for algorithms that +// recursively subdivide cells. +func (rc *RegionCoverer) FastCovering(region Region) CellUnion { + c := rc.newCoverer() + cu := CellUnion(region.CellUnionBound()) + c.normalizeCovering(&cu) + return cu +} + +// normalizeCovering normalizes the "covering" so that it conforms to the current covering +// parameters (MaxCells, minLevel, maxLevel, and levelMod). +// This method makes no attempt to be optimal. In particular, if +// minLevel > 0 or levelMod > 1 then it may return more than the +// desired number of cells even when this isn't necessary. +// +// Note that when the covering parameters have their default values, almost +// all of the code in this function is skipped. +func (c *coverer) normalizeCovering(covering *CellUnion) { + // If any cells are too small, or don't satisfy levelMod, then replace them with ancestors. + if c.maxLevel < maxLevel || c.levelMod > 1 { + for i, ci := range *covering { + level := ci.Level() + newLevel := c.adjustLevel(minInt(level, c.maxLevel)) + if newLevel != level { + (*covering)[i] = ci.Parent(newLevel) + } + } + } + // Sort the cells and simplify them. + covering.Normalize() + + // If there are still too many cells, then repeatedly replace two adjacent + // cells in CellID order by their lowest common ancestor. + for len(*covering) > c.maxCells { + bestIndex := -1 + bestLevel := -1 + for i := 0; i+1 < len(*covering); i++ { + level, ok := (*covering)[i].CommonAncestorLevel((*covering)[i+1]) + if !ok { + continue + } + level = c.adjustLevel(level) + if level > bestLevel { + bestLevel = level + bestIndex = i + } + } + + if bestLevel < c.minLevel { + break + } + (*covering)[bestIndex] = (*covering)[bestIndex].Parent(bestLevel) + covering.Normalize() + } + // Make sure that the covering satisfies minLevel and levelMod, + // possibly at the expense of satisfying MaxCells. + if c.minLevel > 0 || c.levelMod > 1 { + covering.Denormalize(c.minLevel, c.levelMod) + } +} + +// SimpleRegionCovering returns a set of cells at the given level that cover +// the connected region and a starting point on the boundary or inside the +// region. The cells are returned in arbitrary order. +// +// Note that this method is not faster than the regular Covering +// method for most region types, such as Cap or Polygon, and in fact it +// can be much slower when the output consists of a large number of cells. +// Currently it can be faster at generating coverings of long narrow regions +// such as polylines, but this may change in the future. +func SimpleRegionCovering(region Region, start Point, level int) []CellID { + return FloodFillRegionCovering(region, cellIDFromPoint(start).Parent(level)) +} + +// FloodFillRegionCovering returns all edge-connected cells at the same level as +// the given CellID that intersect the given region, in arbitrary order. +func FloodFillRegionCovering(region Region, start CellID) []CellID { + var output []CellID + all := map[CellID]bool{ + start: true, + } + frontier := []CellID{start} + for len(frontier) > 0 { + id := frontier[len(frontier)-1] + frontier = frontier[:len(frontier)-1] + if !region.IntersectsCell(CellFromCellID(id)) { + continue + } + output = append(output, id) + for _, nbr := range id.EdgeNeighbors() { + if !all[nbr] { + all[nbr] = true + frontier = append(frontier, nbr) + } + } + } + + return output +} + +// TODO(roberts): The differences from the C++ version +// finish up FastCovering to match C++ +// IsCanonical +// CanonicalizeCovering +// containsAllChildren +// replaceCellsWithAncestor |