Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

remove geographic_crs and update info model #746

Merged
Merged
Show file tree
Hide file tree
Changes from 1 commit
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
25 changes: 8 additions & 17 deletions rio_tiler/io/base.py
Original file line number Diff line number Diff line change
Expand Up @@ -15,7 +15,7 @@
from rasterio.rio.overview import get_maximum_overview_level
from rasterio.warp import calculate_default_transform, transform_bounds

from rio_tiler.constants import WEB_MERCATOR_TMS, WGS84_CRS
from rio_tiler.constants import WEB_MERCATOR_TMS
from rio_tiler.errors import (
AssetAsBandError,
ExpressionMixingWarning,
Expand All @@ -27,7 +27,7 @@
from rio_tiler.models import BandStatistics, ImageData, Info, PointData
from rio_tiler.tasks import multi_arrays, multi_points, multi_values
from rio_tiler.types import AssetInfo, BBox, Indexes
from rio_tiler.utils import cast_to_sequence, normalize_bounds
from rio_tiler.utils import CRS_to_uri, cast_to_sequence, normalize_bounds


@attr.s
Expand All @@ -48,12 +48,9 @@ class SpatialMixin:
height: Optional[int] = attr.ib(default=None, init=False)
width: Optional[int] = attr.ib(default=None, init=False)

geographic_crs: CRS = attr.ib(init=False, default=WGS84_CRS)

@cached_property
def geographic_bounds(self) -> BBox:
"""Return dataset bounds in geographic_crs."""
if self.crs == self.geographic_crs:
def get_geographic_bounds(self, crs: CRS) -> BBox:
vincentsarago marked this conversation as resolved.
Show resolved Hide resolved
vincentsarago marked this conversation as resolved.
Show resolved Hide resolved
"""Return Geographic Bounds for a Geographic CRS."""
if self.crs == crs:
if self.bounds[1] > self.bounds[3]:
warnings.warn(
"BoundingBox of the dataset is inverted (minLat > maxLat).",
Expand All @@ -69,12 +66,7 @@ def geographic_bounds(self) -> BBox:
return self.bounds

try:
bounds = transform_bounds(
self.crs,
self.geographic_crs,
*self.bounds,
densify_pts=21,
)
bounds = transform_bounds(self.crs, crs, *self.bounds, densify_pts=21)
except: # noqa
warnings.warn(
"Cannot determine bounds in geographic CRS, will default to (-180.0, -90.0, 180.0, 90.0).",
Expand Down Expand Up @@ -1085,9 +1077,8 @@ def _reader(band: str, **kwargs: Any) -> Info:
bands_metadata = multi_values(bands, _reader, **kwargs)

meta = {
"bounds": self.geographic_bounds,
"minzoom": self.minzoom,
"maxzoom": self.maxzoom,
"bounds": self.bounds,
"crs": CRS_to_uri(self.crs) or self.crs.to_wkt(),
}

# We only keep the value for the first band.
Expand Down
16 changes: 8 additions & 8 deletions rio_tiler/io/rasterio.py
Original file line number Diff line number Diff line change
Expand Up @@ -35,7 +35,12 @@
from rio_tiler.io.base import BaseReader
from rio_tiler.models import BandStatistics, ImageData, Info, PointData
from rio_tiler.types import BBox, Indexes, NumType, RIOResampling
from rio_tiler.utils import _validate_shape_input, has_alpha_band, has_mask_band
from rio_tiler.utils import (
CRS_to_uri,
_validate_shape_input,
has_alpha_band,
has_mask_band,
)


@attr.s
Expand All @@ -46,7 +51,6 @@ class Reader(BaseReader):
input (str): dataset path.
dataset (rasterio.io.DatasetReader or rasterio.io.DatasetWriter or rasterio.vrt.WarpedVRT, optional): Rasterio dataset.
tms (morecantile.TileMatrixSet, optional): TileMatrixSet grid definition. Defaults to `WebMercatorQuad`.
geographic_crs (rasterio.crs.CRS, optional): CRS to use as geographic coordinate system. Defaults to WGS84.
colormap (dict, optional): Overwrite internal colormap.
options (dict, optional): Options to forward to low-level reader methods.

Expand Down Expand Up @@ -75,7 +79,6 @@ class Reader(BaseReader):
)

tms: TileMatrixSet = attr.ib(default=WEB_MERCATOR_TMS)
geographic_crs: CRS = attr.ib(default=WGS84_CRS)

colormap: Dict = attr.ib(default=None)

Expand Down Expand Up @@ -177,9 +180,8 @@ def _get_descr(ix):
nodata_type = "None"

meta = {
"bounds": self.geographic_bounds,
"minzoom": self.minzoom,
"maxzoom": self.maxzoom,
"bounds": self.bounds,
"crs": CRS_to_uri(self.crs) or self.crs.to_wkt(),
"band_metadata": [
(f"b{ix}", self.dataset.tags(ix)) for ix in self.dataset.indexes
],
Expand Down Expand Up @@ -596,8 +598,6 @@ class ImageReader(Reader):
tms: TileMatrixSet = attr.ib(init=False)

crs: CRS = attr.ib(init=False, default=None)
geographic_crs: CRS = attr.ib(init=False, default=None)

transform: Affine = attr.ib(init=False)

def __attrs_post_init__(self):
Expand Down
4 changes: 0 additions & 4 deletions rio_tiler/io/stac.py
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,6 @@
from cachetools import LRUCache, cached
from cachetools.keys import hashkey
from morecantile import TileMatrixSet
from rasterio.crs import CRS
from rasterio.transform import array_bounds

from rio_tiler.constants import WEB_MERCATOR_TMS, WGS84_CRS
Expand Down Expand Up @@ -198,7 +197,6 @@ class STACReader(MultiBaseReader):
tms (morecantile.TileMatrixSet, optional): TileMatrixSet grid definition. Defaults to `WebMercatorQuad`.
minzoom (int, optional): Set minzoom for the tiles.
maxzoom (int, optional): Set maxzoom for the tiles.
geographic_crs (rasterio.crs.CRS, optional): CRS to use as geographic coordinate system. Defaults to WGS84.
include_assets (set of string, optional): Only Include specific assets.
exclude_assets (set of string, optional): Exclude specific assets.
include_asset_types (set of string, optional): Only include some assets base on their type.
Expand Down Expand Up @@ -234,8 +232,6 @@ class STACReader(MultiBaseReader):
minzoom: int = attr.ib(default=None)
maxzoom: int = attr.ib(default=None)

geographic_crs: CRS = attr.ib(default=WGS84_CRS)

include_assets: Optional[Set[str]] = attr.ib(default=None)
exclude_assets: Optional[Set[str]] = attr.ib(default=None)

Expand Down
9 changes: 3 additions & 6 deletions rio_tiler/io/xarray.py
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,7 @@
from rio_tiler.io.base import BaseReader
from rio_tiler.models import BandStatistics, ImageData, Info, PointData
from rio_tiler.types import BBox, NoData, WarpResampling
from rio_tiler.utils import _validate_shape_input
from rio_tiler.utils import CRS_to_uri, _validate_shape_input

try:
import xarray
Expand All @@ -43,7 +43,6 @@ class XarrayReader(BaseReader):
Attributes:
dataset (xarray.DataArray): Xarray DataArray dataset.
tms (morecantile.TileMatrixSet, optional): TileMatrixSet grid definition. Defaults to `WebMercatorQuad`.
geographic_crs (rasterio.crs.CRS, optional): CRS to use as geographic coordinate system. Defaults to WGS84.

Examples:
>>> ds = xarray.open_dataset(
Expand All @@ -62,7 +61,6 @@ class XarrayReader(BaseReader):
input: xarray.DataArray = attr.ib()

tms: TileMatrixSet = attr.ib(default=WEB_MERCATOR_TMS)
geographic_crs: CRS = attr.ib(default=WGS84_CRS)

_dims: List = attr.ib(init=False, factory=list)

Expand Down Expand Up @@ -120,9 +118,8 @@ def info(self) -> Info:
metadata = [band.attrs for d in self._dims for band in self.input[d]]

meta = {
"bounds": self.geographic_bounds,
"minzoom": self.minzoom,
"maxzoom": self.maxzoom,
"bounds": self.bounds,
"crs": CRS_to_uri(self.crs),
"band_metadata": [(f"b{ix}", v) for ix, v in enumerate(metadata, 1)],
"band_descriptions": [(f"b{ix}", v) for ix, v in enumerate(bands, 1)],
"dtype": str(self.input.dtype),
Expand Down
10 changes: 2 additions & 8 deletions rio_tiler/models.py
Original file line number Diff line number Diff line change
Expand Up @@ -61,16 +61,10 @@ class Bounds(RioTilerBaseModel):
"""Dataset Bounding box"""

bounds: BoundingBox
crs: str


class SpatialInfo(Bounds):
"""Dataset SpatialInfo"""

minzoom: int
maxzoom: int


class Info(SpatialInfo):
class Info(Bounds):
"""Dataset Info."""

band_metadata: List[Tuple[str, Dict]]
Expand Down
18 changes: 18 additions & 0 deletions rio_tiler/utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -789,3 +789,21 @@ def cast_to_sequence(val: Optional[Any] = None) -> Sequence:
val = (val,)

return val


def CRS_to_uri(crs: CRS) -> Optional[str]:
"""Convert CRS to URI.

Code adapted from https://github.com/developmentseed/morecantile/blob/1829fe12408e4a1feee7493308f3f02257ef4caf/morecantile/models.py#L148-L161
vincentsarago marked this conversation as resolved.
Show resolved Hide resolved
"""
# attempt to grab the authority, version, and code from the CRS
if authority_code := crs.to_authority(confidence_threshold=70):
version = "0"
authority, code = authority_code
# if we have a version number in the authority, split it out
if "_" in authority:
authority, version = authority.split("_")

return f"http://www.opengis.net/def/crs/{authority}/{version}/{code}"

return None
82 changes: 38 additions & 44 deletions tests/test_io_rasterio.py
Original file line number Diff line number Diff line change
@@ -1,18 +1,17 @@
"""tests rio_tiler.io.rasterio.Reader"""

import os
import warnings
from io import BytesIO
from typing import Any, Dict

import attr
import morecantile
import numpy
import pytest
import rasterio
from morecantile import TileMatrixSet
from pyproj import CRS
from rasterio import transform
from rasterio.crs import CRS as rioCRS
from rasterio.features import bounds as featureBounds
from rasterio.io import MemoryFile
from rasterio.vrt import WarpedVRT
Expand Down Expand Up @@ -82,6 +81,10 @@ def test_info_valid():
"""Should work as expected (get file info)"""
with Reader(COG_SCALE) as src:
meta = src.info()
assert meta.bounds == src.bounds == src.dataset.bounds
crs = meta.crs
assert rioCRS.from_user_input(crs) == src.crs

assert meta.scales
assert meta.offsets
assert not meta.colormap
Expand All @@ -92,6 +95,11 @@ def test_info_valid():
assert meta.driver

with Reader(COG_CMAP) as src:
meta = src.info()
assert meta.bounds == src.bounds == src.dataset.bounds
crs = meta.crs
assert rioCRS.from_user_input(crs) == src.crs

assert src.colormap
meta = src.info()
assert meta["colormap"]
Expand All @@ -106,8 +114,7 @@ def test_info_valid():
with Reader(COG_TAGS) as src:
meta = src.info()
assert meta.bounds
assert meta.minzoom
assert meta.maxzoom
assert meta.crs
assert meta.band_descriptions
assert meta.dtype == "int16"
assert meta.colorinterp == ["gray"]
Expand Down Expand Up @@ -894,26 +901,30 @@ def test_nonearthbody():
assert src.maxzoom == 24

# Warns because of zoom level in WebMercator can't be defined
with pytest.warns(UserWarning) as w:
with Reader(COG_EUROPA, geographic_crs=EUROPA_SPHERE) as src:
assert src.info()
assert len(w) == 2
with Reader(COG_EUROPA) as src:
assert src.info()
meta = src.info()
assert meta.bounds == src.bounds == src.dataset.bounds
crs = meta.crs
assert rioCRS.from_user_input(crs) == src.crs

img = src.read()
assert numpy.array_equal(img.data, src.dataset.read(indexes=(1,)))
assert img.width == src.dataset.width
assert img.height == src.dataset.height
assert img.count == src.dataset.count
assert src.get_geographic_bounds(EUROPA_SPHERE)

img = src.read()
assert numpy.array_equal(img.data, src.dataset.read(indexes=(1,)))
assert img.width == src.dataset.width
assert img.height == src.dataset.height
assert img.count == src.dataset.count

img = src.preview()
assert img.bounds == src.bounds
img = src.preview()
assert img.bounds == src.bounds

part = src.part(src.bounds, bounds_crs=src.crs)
assert part.bounds == src.bounds
part = src.part(src.bounds, bounds_crs=src.crs)
assert part.bounds == src.bounds

lon = (src.bounds[0] + src.bounds[2]) / 2
lat = (src.bounds[1] + src.bounds[3]) / 2
assert src.point(lon, lat, coord_crs=src.crs).data[0] is not None
lon = (src.bounds[0] + src.bounds[2]) / 2
lat = (src.bounds[1] + src.bounds[3]) / 2
assert src.point(lon, lat, coord_crs=src.crs).data[0] is not None

with pytest.warns(UserWarning):
europa_crs = CRS.from_authority("ESRI", 104915)
Expand All @@ -923,7 +934,7 @@ def test_nonearthbody():
matrix_scale=[2, 1],
)

with Reader(COG_EUROPA, tms=tms, geographic_crs=EUROPA_SPHERE) as src:
with Reader(COG_EUROPA, tms=tms) as src:
assert src.info()
assert src.minzoom == 4
assert src.maxzoom == 6
Expand Down Expand Up @@ -958,28 +969,8 @@ def test_nonearth_custom():
geographic_crs=MARS2000_SPHERE,
)

@attr.s
class MarsReader(Reader):
"""Use custom geographic CRS."""

geographic_crs: rasterio.crs.CRS = attr.ib(
init=False,
default=rasterio.crs.CRS.from_proj4("+proj=longlat +R=3396190 +no_defs"),
)

with warnings.catch_warnings():
with MarsReader(COG_MARS, tms=mars_tms) as src:
assert src.geographic_bounds[0] > -180

with warnings.catch_warnings():
with Reader(
COG_MARS,
tms=mars_tms,
geographic_crs=rasterio.crs.CRS.from_proj4(
"+proj=longlat +R=3396190 +no_defs"
),
) as src:
assert src.geographic_bounds[0] > -180
with Reader(COG_MARS, tms=mars_tms) as src:
assert src.get_geographic_bounds(MARS2000_SPHERE)[0] > -180


def test_tms_tilesize_and_zoom():
Expand Down Expand Up @@ -1121,7 +1112,10 @@ def test_inverted_latitude():
"""Test working with inverted Latitude."""
with pytest.warns(UserWarning):
with Reader(COG_INVERTED) as src:
assert src.geographic_bounds[1] < src.geographic_bounds[3]
assert (
src.get_geographic_bounds(WGS84_CRS)[1]
< src.get_geographic_bounds(WGS84_CRS)[3]
)

with pytest.warns(UserWarning):
with Reader(COG_INVERTED) as src:
Expand Down
12 changes: 7 additions & 5 deletions tests/test_io_xarray.py
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,7 @@
import pytest
import rioxarray
import xarray
from rasterio.crs import CRS as rioCRS

from rio_tiler.errors import InvalidGeographicBounds, MissingCRS
from rio_tiler.io import XarrayReader
Expand All @@ -28,9 +29,11 @@ def test_xarray_reader():

data.rio.write_crs("epsg:4326", inplace=True)
with XarrayReader(data) as dst:
assert dst.minzoom == dst.maxzoom == 0
info = dst.info()
assert info.minzoom == 0
assert info.maxzoom == 0
assert info.bounds == dst.bounds
crs = info.crs
assert rioCRS.from_user_input(crs) == dst.crs
assert info.band_metadata == [("b1", {})]
assert info.band_descriptions == [("b1", "2022-01-01T00:00:00.000000000")]
assert info.height == 33
Expand Down Expand Up @@ -123,9 +126,8 @@ def test_xarray_reader():

data.rio.write_crs("epsg:4326", inplace=True)
with XarrayReader(data) as dst:
info = dst.info()
assert info.minzoom == 5
assert info.maxzoom == 7
assert dst.minzoom == 5
assert dst.maxzoom == 7


def test_xarray_reader_external_nodata():
Expand Down
Loading