Skip to content

The ARD Grid

Maxar ARD uses a grid system based on UTM zones. Each zone is divided using a quadkey scheme that divides the zone into useful grid cells.

At the 12th zoom level each cell is 5 kilometers square. We use these cells as a boundary to subset image strips into tiles. The image tiles are Cloud Optimized Geotiffs stored in S3. The tiffs use a normalized resolution of approximately 30cm. The size and resolution of the tiffs maps to 16,384 pixels inside of each cell (the COGS are generated with an additional 512 pixel apron outside of the cell boundary). Because 16,384 is a power of two, it can be further subdivided with quadkeys along pixel boundaries.

At the 18th zoom level we have a cell that is 256 pixels in size. This is a convenient size for a chip. Each chip therefore has a unique quadkey describing its bounding box. Fetching image data for this cell will alway retrieve the same square on the earth. Because 256 is also a power of two, that means every pixel also has a unique quadkey identifier on the globe.

maxar-ard-grid

maxar-ard-grid is a Python package for working with the ARD grid. It is included as a dependency in the max-ard SDK but can also be used standalone.

Note: Prior to July 15th 2020 this package was called maxar-canvas-grid. The name has been changed and the package is now available in PyPI. maxar-canvas-grid will remain available at the 1.1.0 version from packages.ard.maxar.com until January 1, 2023. It will not receive any updates however, so we recommend anyone using the old package update their dependencies and imports.

Installation

Install using pip:

pip install maxar-ard-grid

The Basics of the ARD Grid

Finding grid cells

To find out which grid cells cover an area, use the covers() function:

>>> from maxar_ard_grid import covers
>>> geom = 'POINT(-120 40)'
>>> covers(geom)

<generator object _covers_geom at 0x107ef6890>

The covers function returns a generator of Cell objects. Typically one would would iterate over the generator:

for cell in covers(geom):
    print(cell)

<Cell Z10-120020112013>

To turn the iterator into a list, just call list() on it:

cells = list(covers(geom))
print(len(cells))

1

Unless specified the Cells will be generated at grid zoom level 12 to match the ARD tiles. Passing an optional zoom parameter will return cells at different zoom levels. Therefore these two are equivalent:

>>> covers(geom)
>>> covers(geom, zoom=12)

The following are acceptable formats for the input geometry:

  • Objects using the Python Geospatial Interface
  • Shapely geometries
  • Cell IDs
  • Cell Objects
  • GeoJSON-like geometry dicts
  • GeoJSON-like feature dicts
  • Well Know Text (WKT)

Initalizing a Cell with a Cell returns the same cell.

Geometries are assumed to be in the WGS84 coordinate system. If the geometry is in a different coordinate system or projection, it can be specified using the src_proj keyword. Acceptable inputs are PyProj CRS objects or EPSG codes (as strings or integers):

geom = 'POINT(243900 4432069)'
list(covers(geom, src_proj=32611)) # 32611 is UTM 11N

[<Cell Z10-120020112013>]

Working with Cells

Each cell in a zone's grid has a unique quadkey identifier. When combined with the the zone's UTM number the result is a unique cell on Earth. These cell IDs use the following format:

Z{zone}-{quadkey}

An example is Z08-03201333302.

Note that the zone number is zero padded to two digits. Quadkeys can start with the digit 0 so they should always be represented as strings. The length of a quadkey is equal to the cell's zoom level.

A specific cell can be initialized using either a Cell ID, or a combination of quadkey and zone:

from maxar_ard_grid import Cell
c = Cell('Z08-032013333021')
c = Cell('032013333021', zone=8)

Identifying attributes:

c.id

'Z08-032013333021'

c.quadkey

'032013333021'

c.zone

8

c.zoom

12

Geospatial attributes

Cell objects have several attributes and methods for geometric and geospatial functionality:

c.proj # Projection EPSG code

'EPSG:32608'

c.bounds # bounds in UTM

[-3375000.0, 1945000.0, -3370000.0, 1950000.0]

# Python Geospatial Interface support
from shapely.geometry import shape
shape(c) # a Polygon in UTM


c.geom_WGS84 # Polygon in WGS84 coordinates


c.to_feature() # GeoJSON-like Python dictionary

{'type': 'Feature',
 'properties': {'id': 'Z08-032013333021'},
 'geometry': {'type': 'Polygon',
  'coordinates': (((-169.06140140640161, 14.697247231623834),
    (-169.06801910943005, 14.734648088759604),
    (-169.10643859320024, 14.72819673954847),
    (-169.09981733357938, 14.690812971571576),
    (-169.06140140640161, 14.697247231623834)),)}}

c.to_GeoJSON() # GeoJSON feature string

'{"type": "Feature", "properties": {"id": "Z08-032013333021"}, "geometry": {"type": "Polygon", "coordinates": [[[-169.06140140640161, 14.697247231623834], [-169.06801910943005, 14.734648088759604], [-169.10643859320024, 14.72819673954847], [-169.09981733357938, 14.690812971571576], [-169.06140140640161, 14.697247231623834]]]}}'

c.parent

<Cell Z08-03201333302>

c.get_parent_at_zoom(9)

<Cell Z08-032013333>

c.children

[<Cell Z08-0320133330210>,
 <Cell Z08-0320133330211>,
 <Cell Z08-0320133330212>,
 <Cell Z08-0320133330213>]

c.get_children_at_zoom(18) # returns a generator of cells 

<generator object Cell.get_children_at_zoom at 0x107ef6820>

c.neighbor('E') # cell to the east

<Cell Z08-032013333030>

c.neighbors

{'N': <Cell Z08-032013333003>,
 'NE': <Cell Z08-032013333012>,
 'E': <Cell Z08-032013333030>,
 'SE': <Cell Z08-032013333032>,
 'S': <Cell Z08-032013333023>,
 'SW': <Cell Z08-032013333022>,
 'W': <Cell Z08-032013333020>,
 'NW': <Cell Z08-032013333002>}

Overloaded operators include == and in. Cell objects are equivalent if they represent the same cell on earth. in overloads __contains__ to check if a cell is a descendent of a given cell, which is also equivalent to being spatially contained.

Cell('Z08-03201333302') == Cell('03201333302', zone=8)

True

Cell('Z08-03201333302') in Cell('Z08-032013')

True

Cautions and Caveats

It is possible to get confusing results when generating covering cells. If you get more or fewer cells than you expect, it could be due to one of these "gotchas":

Lower dimension intersections.

A geometry technically intersects a cell if it shares an edge or even a vertex. This library tries to be smart by checking intersections using an ever-so-slightly smaller boundary so that the cells covered by a cell are only said cell, instead of the cell and its 8 neighbors that share edges and vertices.

If you find yourself getting extra cells, you may need to verify that your geometry overlaps the cell by checking if the intersection of the two geometries have a non-zero area.

The Earth isn't flat.

Cells are only represented as 4 vertices joined by straight lines in the UTM projection. In many software implementations, reprojection is only carried out by reprojecting the vertex points. Their connecting edge however is still represented as a straight line! If you are needing the edges of a cell to represent the curve properly in geographic coordinate systems you will need to densify the edges to add intermediate vertices. When reprojected they will better represent the curvature of the shape.

This is especially noticable in trying to use the shape of a higher zoom level cell to select lower zoom level cells. The straight edge spanning the longer distance does not track the individual edges of the smaller cells.

Back to top