rio-tiler 2.0b7

Rasterio plugin to read mercator tiles from Cloud Optimized GeoTIFF.

Rio-tiler

Rasterio plugin to read mercator tiles from Cloud Optimized GeoTIFF.

Install

You can install rio-tiler using pip

$ pip install -U pip
$ pip install rio-tiler --pre # version 2.0 is in development

or install from source:

$ git clone https://github.com/cogeotiff/rio-tiler.git
$ cd rio-tiler
$ pip install -U pip
$ pip install -e .

Usage

The rio_tiler module can create mercator tiles from any raster source supported by Rasterio/GDAL (i.e. local files, http, s3, gcs etc.). Additional method are availables (see COGReader)

Read a tile from a file

from rio_tiler.io import COGReader

with COGReader("http://oin-hotosm.s3.amazonaws.com/5a95f32c2553e6000ce5ad2e/0/10edab38-1bdd-4c06-b83d-6e10ac532b7d.tif") as cog:
    tile, mask = cog.tile(691559, 956905, 21, tilesize=256)

print(tile.shape)
> (3, 256, 256)

print(mask.shape)
> (256, 256)

Render the array as an image (PNG/JPEG)

from rio_tiler.utils import render

buffer = render(tile, mask=mask) # this returns a buffer (PNG by default)

Rescale non-byte data and/or apply colormap

from rio_tiler.colormap import cmap
from rio_tiler.utils import linear_rescale

# Rescale the tile array only where mask is valid and cast it to byte
tile = numpy.where(
    mask,
    linear_rescale(tile, in_range=(0, 1000), out_range=[0, 255]),
    0
).astype(numpy.uint8)

cm = cmap.get("viridis")

buffer = render(tile, mask=mask, colormap=cm)

Use creation options to match mapnik defaults.

from rio_tiler.utils import render
from rio_tiler.profiles import img_profiles

options = img_profiles.get("webp")
buffer = render(tile, mask=mask, img_format="webp", **options)

Write image to file

with open("my.png", "wb") as f:
  f.write(buffer)

COGReader

class COGReader:
    """
    Cloud Optimized GeoTIFF Reader.

    Examples
    --------
    with COGReader(src_path) as cog:
        cog.tile(...)

    # Set global options
    with COGReader(src_path, unscale=True, nodata=0) as cog:
        cog.tile(...)

    with rasterio.open(src_path) as src_dst:
        with WarpedVRT(src_dst, ...) as vrt_dst:
            with COGReader(None, dataset=vrt_dst) as cog:
                cog.tile(...)

    with rasterio.open(src_path) as src_dst:
        with COGReader(None, dataset=src_dst) as cog:
            cog.tile(...)

    Attributes
    ----------
    filepath: str
        Cloud Optimized GeoTIFF path.
    dataset: rasterio.DatasetReader, optional
        Rasterio dataset.

    Properties
    ----------
    minzoom: int
        COG minimum zoom level.
    maxzoom: int
        COG maximum zoom level.
    bounds: tuple[float]
        COG bounds in WGS84 crs.
    center: tuple[float, float, int]
        COG center + minzoom
    colormap: dict
        COG internal colormap.

    Methods
    -------
    tile(0, 0, 0, indexes=(1,2,3), expression="
B1/B2", tilesize=512, resampling_methods="nearest")
        Read a map tile from the COG.
    part((0,10,0,10), indexes=(1,2,3,), expression="
B1/B20", max_size=1024)
        Read part of the COG.
    preview(max_size=1024)
        Read preview of the COG.
    point((10, 10), indexes=1)
        Read a point value from the COG.
    info: dict
        General information about the COG (datatype, indexes, ...)
    stats(pmin=5, pmax=95)
        Get Raster statistics.
    meta(pmin=5, pmax=95)
        Get info + raster statistics
    """

Properties

• dataset: Return the rasterio dataset
• colormap: Return the dataset's internal colormap
• minzoom: Return minimum Mercator Zoom
• maxzoom: Return maximum Mercator Zoom
• bounds: Return the dataset bounds in WGS84
• center: Return the center of the dataset + minzoom
• spatial_info: Return the bounds, center and zoom infos

Methods

• tile(): Read map tile from a raster

with COGReader("myfile.tif") as cog:
    tile, mask = cog.tile(1, 2, 3, tilesize=256)

# With indexes
with COGReader("myfile.tif") as cog:
    tile, mask = cog.tile(1, 2, 3, tilesize=256, indexes=1)

# With expression
with COGReader("myfile.tif"s) as cog:
    tile, mask = cog.tile(1, 2, 3, tilesize=256, expression="B1/B2")

• part(): Read part of a raster

with COGReader("myfile.tif") as cog:
    data, mask = cog.part((10, 10, 20, 20))

# Limit output size (default is set to 1024)
with COGReader("myfile.tif") as cog:
    data, mask = cog.part((10, 10, 20, 20), max_size=2000)

# Read high resolution
with COGReader("myfile.tif") as cog:
    data, mask = cog.part((10, 10, 20, 20), max_size=None)

# With indexes
with COGReader("myfile.tif") as cog:
     data, mask = cog.part((10, 10, 20, 20), indexes=1)

# With expression
with COGReader("myfile.tif") as cog:
    data, mask = cog.part((10, 10, 20, 20), expression="B1/B2")

• preview(): Read a preview of a raster

with COGReader("myfile.tif") as cog: 
    data, mask = cog.preview()

# With indexes
with COGReader("myfile.tif") as cog: 
    data, mask = cog.preview(indexes=1)

# With expression
with COGReader("myfile.tif") as cog: 
    data, mask = cog.preview(expression="B1+2,B1*4")

• point(): Read point value of a raster

with COGReader("myfile.tif") as cog: 
    print(cog.point(-100, 25))

# With indexes
with COGReader("myfile.tif") as cog: 
    print(cog.point(-100, 25, indexes=1)) 
[1]

# With expression
with COGReader("myfile.tif") as cog: 
    print(cog.point(-100, 25, expression="B1+2,B1*4"))
[3, 4]

• info(): Return simple metadata about the dataset

with COGReader("myfile.tif") as cog:
    print(cog.info())
{
    "bounds": [-119.05915661478785, 13.102845359730287, -84.91821332299578, 33.995073647795806],
    "center": [-101.98868496889182, 23.548959503763047, 3],
    "minzoom": 3,
    "maxzoom": 12,
    "band_metadata": [[1, {}]],
    "band_descriptions": [[1,"band1"]],
    "dtype": "int8",
    "colorinterp": ["palette"],
    "nodata_type": "Nodata",
    "colormap": {
        "0": [0, 0, 0, 0],
        "1": [0, 61, 0, 255],
        ...
    }
}

• stats(): Return image statistics (Min/Max/Stdev)

with COGReader("myfile.tif") as cog:
    print(cog.stats())
{
    "1": {
        "pc": [1, 16],
        "min": 1,
        "max": 18,
        "std": 4.069636227214257,
        "histogram": [
            [...],
            [...]
        ]
    }
}

• metadata(): Return COG info + statistics

with COGReader("myfile.tif") as cog:
    print(cog.metadata())
{
    "bounds": [-119.05915661478785, 13.102845359730287, -84.91821332299578, 33.995073647795806],
    "center": [-101.98868496889182, 23.548959503763047, 3],
    "minzoom": 3,
    "maxzoom": 12,
    "band_metadata": [[1, {}]],
    "band_descriptions": [[1,"band1"]],
    "dtype": "int8",
    "colorinterp": ["palette"],
    "nodata_type": "Nodata",
    "colormap": {
        "0": [0, 0, 0, 0],
        "1": [0, 61, 0, 255],
        ...
    }
    "statistics" : {
        1: {
            "pc": [1, 16],
            "min": 1,
            "max": 18,
            "std": 4.069636227214257,
            "histogram": [
                [...],
                [...]
            ]
        }
    }
}

Global Options

COGReader accept several options which will be forwarded to the rio_tiler.reader._read function (low level function accessing the data):

• nodata: Overwrite the nodata value (or set if not present)
• unscale: Apply internal rescaling factors
• vrt_options: Pass WarpedVRT Option (see: https://gdal.org/api/gdalwarp_cpp.html?highlight=vrt#_CPPv415GDALWarpOptions)

Note: Those options could already be passed on each method call.

with COGReader("my_cog.tif", nodata=0) as cog:
    tile, mask = cog.tile(1, 1, 1)

# is equivalent to 

with COGReader("my_cog.tif") as cog:
    tile, mask = cog.tile(1, 1, 1, nodata=0)

STACReader

In rio-tiler v2, we added a rio_tiler.io.STACReader to allow tile/metadata fetching of assets withing a STAC item. The STACReader objects has the same properties/methods as the COGReader.

from typing import Dict 
from rio_tiler.io import STACReader

with STACReader(
    "https://1tqdbvsut9.execute-api.us-west-2.amazonaws.com/v0/collections/sentinel-s2-l2a-cogs/items/S2A_34SGA_20200318_0_L2A", 
    exclude_assets={"thumbnail"} 
) as stac: 
    print(stac.bounds)
    print(stac.assets)

> [23.293255090449595, 31.505183020453355, 24.296453548295318, 32.51147809805106]
> ['overview', 'visual', 'B01', 'B02', 'B03', 'B04', 'B05', 'B06', 'B07', 'B08', 'B8A', 'B09', 'B11', 'B12', 'AOT', 'WVP', 'SCL']

# Name of assets to read
assets = ["B01", "B02"]

with STACReader(
    "https://1tqdbvsut9.execute-api.us-west-2.amazonaws.com/v0/collections/sentinel-s2-l2a-cogs/items/S2A_34SGA_20200318_0_L2A", 
    exclude_assets={"thumbnail"} 
) as stac:
    tile, mask = stac.tile(145, 103, 8, tilesize=256, assets=assets)

print(tile.shape)
> (2, 256, 256)

# With expression
with STACReader(
    "https://1tqdbvsut9.execute-api.us-west-2.amazonaws.com/v0/collections/sentinel-s2-l2a-cogs/items/S2A_34SGA_20200318_0_L2A", 
    exclude_assets={"thumbnail"} 
) as stac:
    tile, mask = stac.tile(145, 103, 8, tilesize=256, expression="B01/B02")

print(tile.shape)
> (1, 256, 256)

Note: STACReader is based on rio_tiler.io.base.MultiBaseReader class.

Working with multiple assets

Mosaic

Starting in rio-tiler 2.0, we've transfered the rio-tiler-mosaic plugin to be a rio-tiler submodule.

from rio_tiler.io import COGReader
from rio_tiler.mosaic import mosaic_reader
from rio_tiler.mosaic.methods import defaults


def tiler(src_path: str, *args, **kwargs) -> Tuple[numpy.ndarray, numpy.ndarray]:
    with COGReader(src_path) as cog:
        return cog.tile(*args, **kwargs)

assets = ["mytif1.tif", "mytif2.tif", "mytif3.tif"]
(tile, mask), assets_used = mosaic_reader(assets, tiler, 1, 1, 1)

Learn more about rio_tiler.mosaic in doc/mosaic.md.

Notebook: WorkingWithMosaic

Merge assets

rio_tiler.io.cogeo submodule has multi_* functions (tile, part, preview, point, metadata, info, stats) allowing to fetch and merge info/data from multiple dataset (think about multiple bands stored in separated files).

from typing import Dict 
from rio_tiler.io.cogeo import multi_tile

assets = ["b1.tif", "b2.tif", "b3.tif"]
tile, mask = multi_tile(assets, x, y, z, tilesize=256)

print(tile.shape)
> (3, 256, 256)

# Others
metadata = multi_info(assets)
stats = multi_stats(assets, pmin=2, pmax=98, ...)
metadata = multi_metadata(assets, pmin=2, pmax=98, ...)
values = multi_points(assets, lon, lat, ...)
data, mask = multi_part(assets, bbox, ...)
data, mask = multi_preview(assets, ...)

You can also use rio_tiler.io.base.MultiBaseReader to build a custom asset reader:

import attr
from rio_tiler.io.base import MultiBaseReader
from rio_tiler.io import COGReader, BaseReader


# CustomReader is a subclass of MultiBaseReader.
# To ease the creation of the class and because MultiBaseReader is built with `attr`
# we also need to add the `@attr.s` wrapper on top of our custom class.
@attr.s
class CustomReader(MultiBaseReader):

    directory: str = attr.ib() # required arg
    reader: Type[BaseReader] = attr.ib(default=COGReader) # the default reader is COGReader

    def __enter__(self):
        # List files in directory
        dirs = os.listdir(self.directory) 

        # get list of tifs
        tiff = [f for f in dirs if f.endswith(".tif")]

        # create list of assets names - REQUIRED
        self.assets = [os.path.basename(f).split(".")[0] for f in tiff]
        
        # `self.bounds` needs to be set! - REQUIRED
        with self.reader(tiff[0]) as cog:
            self.bounds = cog.bounds

        return self

    def _get_asset_url(self, asset: str) -> str:
        """Validate asset names and return asset's url."""
        if asset not in self.assets:
            raise InvalidAssetName(f"{asset} is not valid")
        
        return os.path.join(self.directory, f"{asset}.tif")

# we have a directoty with "b1.tif", "b2.tif", "b3.tif"
with CustomReader("my_dir/") as cr:
    print(cr.assets)
    tile, mask = cr.tile(x, y, z, assets="b1")

> ["b1", "b2", "b3"]

print(tile.shape)
> (3, 256, 256)

Reading asset with a GeoJSON Polygon

Natively rio-tiler support mostly bbox reading. Using GDALWarpVRT Cutline option, it's possible to read a dataset for a given polygon.

from rio_tiler.io import COGReader
from rio_tiler.utils import create_cutline
from rasterio.features import bounds as featureBounds

feat =     {
    "type": "Feature",
    "properties": {},
    "geometry": {
    "type": "Polygon",
    "coordinates": [
        [
        [-52.6025390625, 73.86761239709705],
        [-52.6025390625, 73.59679245247814],
        [-51.591796875, 73.60299628304274],
        [-51.591796875, 73.90420357134279],
        [-52.4267578125, 74.0437225981325],
        [-52.6025390625, 73.86761239709705]
        ]
    ]
    }
}

# Get BBOX of the polygon
bbox = featureBounds(feat)

# Use COGReader to open and read the dataset
with COGReader("my_tif.tif") as cog:
    # Create WTT Cutline
    cutline = create_cutline(cog.dataset, feat, geometry_crs="epsg:4326")

    # Read part of the data (bbox) and use the cutline to mask the data
    data, mask = cog.part(bbox, vrt_options={'cutline': cutline})

The previous example uses the .part method but any method that uses the rio_tiler.reader._read function will accept the cutline options.

bbox = featureBounds(feat)

# Use COGReader to open and read the dataset
with COGReader("my_tif.tif") as cog:
    # Create WTT Cutline
    cutline = create_cutline(cog.dataset, feat, geometry_crs="epsg:4326")

    # Get a preview of the whole geotiff but use the cutline to mask the data
    data, mask = cog.preview(vrt_options={'cutline': cutline})

    # Read a mercator tile and use the cutline to mask the data
    data, mask = cog.tile(1, 1, 1, vrt_options={'cutline': cutline})

    # Get image statistics over a bbox and use the cutline as mask
    stats = cog.stats(bounds=bbox, vrt_options={'cutline': cutline})

Partial reading on Cloud hosted dataset

Rio-tiler perform partial reading on local or distant dataset, which is why it will perform best on Cloud Optimized GeoTIFF (COG). It's important to note that Sentinel-2 scenes hosted on AWS are not in Cloud Optimized format but in JPEG2000. When performing partial reading of JPEG2000 dataset GDAL (rasterio backend library) will need to make a lot of GET requests and transfer a lot of data.

Ref: Do you really want people using your data blog post.

Create an AWS Lambda package

The easiest way to make sure the package will work on AWS is to use docker

FROM lambci/lambda:build-python3.7

ENV LANG=en_US.UTF-8 LC_ALL=en_US.UTF-8 CFLAGS="--std=c99"

RUN pip3 install rio-tiler --no-binary numpy -t /tmp/python -U

RUN cd /tmp/python && zip -r9q /tmp/package.zip *

Ref: https://github.com/vincentsarago/simple-rio-lambda

Mission Specific tiler

In rio-tiler v2 we choosed to remove the mission specific tilers (Sentinel2, Sentinel1, Landsat8 and CBERS). Those are now in a specific plugin: rio-tiler-pds.

Plugins

• rio-tiler-mvt: Create Mapbox Vector Tile from numpy array (tile/mask)
• rio-tiler-crs: Create Map Tiles using other TileMatrixSets
• rio-viz: Visualize Cloud Optimized GeoTIFF in browser locally

Implementations

• CosmiQ/solaris
• cogeo-tiler
• titiler

Contribution & Development

Issues and pull requests are more than welcome.

dev install

$ git clone https://github.com/cogeotiff/rio-tiler.git
$ cd rio-tiler
$ pip install -e .[dev]

Python3.7 only

This repo is set to use pre-commit to run isort, flake8, pydocstring, black ("uncompromising Python code formatter") and mypy when committing new code.

$ pre-commit install

License

See LICENSE.txt

Authors

The rio-tiler project was begun at Mapbox and has been transferred in January 2019.

See AUTHORS.txt for a listing of individual contributors.

Changes

See CHANGES.txt.