Using Python S2/S2sphere library - find all s2 cells of a particular level with in a circle(lat, long and radius is given)

Question

What S2Region and how should I use to get all s2 cells at certain parent level(let say 9) covered by the circle drawn from given lat, long and radius. Below is an example which use python s2 library for getting all cells under a rectangle.

region_rect = S2LatLngRect(
              S2LatLng.FromDegrees(-51.264871, -30.241701),
              S2LatLng.FromDegrees(-51.04618, -30.000003))
coverer = S2RegionCoverer()
coverer.set_min_level(8)
coverer.set_max_level(15)
coverer.set_max_cells(500)
covering = coverer.GetCovering(region_rect)

source of example http://blog.christianperone.com/2015/08/googles-s2-geometry-on-the-sphere-cells-and-hilbert-curve/

I am looking for something like

region_circle = S2latLangCircle(lat,lang,radius)

I find answer of this question for google s2 library implemented in c++ Using google s2 library - find all s2 cells of a certain level within the circle, given lat/lng and radius in miles/km but I need this in python.

Thanks

Marna

Answer 1

With the help of link, I worked out for python solution.

I am using python s2sphere library.

earthCircumferenceMeters = 1000 * 40075.017
def earthMetersToRadians(meters):
    return (2 * math.pi) * (float(meters) / 
    const.earthCircumferenceMeters)


def getCoveringRect(lat, lng, radius, parent_level):
    radius_radians = earthMetersToRadians(radius)
    latlng = LatLng.from_degrees(float(lat), 
             float(lng)).normalized().to_point()
    region = Cap.from_axis_height(latlng, 
    (radius_radians*radius_radians)/2)
    coverer = RegionCoverer()
    coverer.min_level = int(parent_level)
    coverer.max_level = int(parent_level)
    coverer.max_cells = const.MAX_S2_CELLS
    covering = coverer.get_covering(region)
    s2_rect = []
    for cell_id in covering:
    new_cell = Cell(cell_id)
    vertices = []
    for i in range(4):
        vertex = new_cell.get_vertex(i)
        latlng = LatLng.from_point(vertex)
        vertices.append((math.degrees(latlng.lat().radians),
                         math.degrees(latlng.lng().radians)))
    s2_rect.append(vertices)
    return s2_rect

getCoveringRect method returns all s2 cells(Rectangle boundary) at given parent level which is covered by circle drawn from given lat, long as center and given radius

Marna

Answer 2

Here is a Go example how to get covering cells

import (
    "math"
    "sort"
    "strconv"

    "github.com/golang/geo/s2"
)

const (
    earthRadiusInMeter = 1000 * 6371.393 // earth radius is 6371km
)

// S=4πR²，s2 regards surface is 4π，that is R=1
func getS2EarthSurfaceArea(radius float64) float64 {
    area := math.Pi * radius * radius / (earthRadiusInMeter * earthRadiusInMeter)
    return area
}


func GetCellIDs(lng, lat, radius float64) []string {
    point := s2.PointFromLatLng(s2.LatLngFromDegrees(lat, lng))
    area := getS2EarthSurfaceArea(radius)
    _cap := s2.CapFromCenterArea(point, area)
    _cover := s2.RegionCoverer{
        MinLevel: 13,
        MaxLevel: 13,
        LevelMod: 1,
        MaxCells: 16,
    }
    cellUnion := _cover.Covering(_cap)
    stringCellIDs := make([]string, 0, len(cellUnion))
    for _, c := range cellUnion {
        stringCellIDs = append(stringCellIDs, c.ToToken())
    }
    return stringCellIDs
}

Gastongastralgia

Answer 3

I'm not sure if the formula used by Guarav is right.

First, function earthMetersToRadians does not return radians, it just computes (2*pi*r) / (2*pi*R) = r/R where R denotes the earth' radius. From that, it computes height = (r/R)^2/2, and I'm not sure where this formula comes from.

From the formulae of a spherical cap, we have height = 1 - cos(theta) where theta = arcsin(r/R) in our case. Together have height = 1 - cos(arcsin(r/R)) which can be computed as height = 1 - sqrt(1 - (r/R)^2).

Note, however, that both formulas are very close, so in practical cases they are pretty much the same, especially if you run an S2Coverer on your cap afterwards.

Tackling