Skip to content

Commit

Permalink
OGC Coverages: add CoverageScaling class
Browse files Browse the repository at this point in the history
This is a helper class for implementing coverage scaling operations.
Not yet fully implemented and not yet in use.
  • Loading branch information
pont-us committed Dec 8, 2023
1 parent 3204918 commit d47f59e
Show file tree
Hide file tree
Showing 3 changed files with 192 additions and 9 deletions.
56 changes: 56 additions & 0 deletions test/webapi/ows/coverages/test_scaling.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,56 @@
import unittest

import pyproj

import xcube.core.new
from xcube.server.api import ApiError
from xcube.webapi.ows.coverages.request import CoveragesRequest
from xcube.webapi.ows.coverages.scaling import CoverageScaling


class ScalingTest(unittest.TestCase):

def setUp(self):
self.epsg4326 = pyproj.CRS('EPSG:4326')
self.ds = xcube.core.new.new_cube()

def test_default_scaling(self):
scaling = CoverageScaling(CoveragesRequest({}), self.epsg4326, self.ds)
self.assertEqual((1, 1), scaling.scale)

def test_no_data(self):
with self.assertRaises(ApiError.NotFound):
CoverageScaling(
CoveragesRequest({}), self.epsg4326,
self.ds.isel(lat=slice(0, 0))
)

def test_crs_fallbacks(self):
# TODO: implement me
pass

def test_scale_factor(self):
scaling = CoverageScaling(
CoveragesRequest({'scale-factor': ['2']}),
self.epsg4326,
self.ds
)
self.assertEqual((2, 2), scaling.scale)

def test_scale_axes(self):
scaling = CoverageScaling(
CoveragesRequest({'scale-axes': ['Lat(3),Lon(1.2)']}),
self.epsg4326,
self.ds
)
self.assertEqual((1.2, 3), scaling.scale)
self.assertEqual((300, 60), scaling.size)

def test_scale_size(self):
scaling = CoverageScaling(
CoveragesRequest({'scale-size': ['Lat(90),Lon(240)']}),
self.epsg4326,
self.ds
)
self.assertEqual((240, 90), scaling.size)
self.assertEqual((1.5, 2), scaling.scale)
23 changes: 14 additions & 9 deletions xcube/webapi/ows/coverages/controllers.py
Original file line number Diff line number Diff line change
Expand Up @@ -206,8 +206,8 @@ def get_coverage_data(

def _assert_coverage_size_ok(ds: xr.Dataset, scale_factor: float):
size_limit = 4000 * 4000 # TODO make this configurable
h_dim = _get_h_dim(ds)
v_dim = _get_v_dim(ds)
h_dim = get_h_dim(ds)
v_dim = get_v_dim(ds)
for d in h_dim, v_dim:
size = ds.dims[d]
if size == 0:
Expand Down Expand Up @@ -338,8 +338,8 @@ def _apply_geographic_subsetting(

# 6. Apply the dataset-native bbox using sel, making sure that y/latitude
# slice has the same ordering as the corresponding co-ordinate.
h_dim = _get_h_dim(ds)
v_dim = _get_v_dim(ds)
h_dim = get_h_dim(ds)
v_dim = get_v_dim(ds)
ds = ds.sel(
indexers={
h_dim: slice(bbox_native_crs[0], bbox_native_crs[2]),
Expand All @@ -354,12 +354,17 @@ def _apply_geographic_subsetting(


def get_bbox_from_ds(ds: xr.Dataset):
h, v = ds[_get_h_dim(ds)], ds[_get_v_dim(ds)]
h, v = ds[get_h_dim(ds)], ds[get_v_dim(ds)]
bbox = list(map(float, [h[0], v[0], h[-1], v[-1]]))
_ensure_bbox_y_ascending(bbox)
return bbox


def apply_scaling(gm: GridMapping, ds: xr.Dataset, request: CoveragesRequest):
# TODO: implement me
pass

Check warning on line 365 in xcube/webapi/ows/coverages/controllers.py

View check run for this annotation

Codecov / codecov/patch

xcube/webapi/ows/coverages/controllers.py#L365

Added line #L365 was not covered by tests


def _find_geographic_parameters(
names: list[str],
) -> tuple[Optional[str], Optional[str]]:
Expand Down Expand Up @@ -396,8 +401,8 @@ def _apply_bbox(
)
_ensure_bbox_y_ascending(bbox, always_xy or is_xy_order(bbox_crs))
bbox = transformer.transform_bounds(*bbox)
h_dim = _get_h_dim(ds)
v_dim = _get_v_dim(ds)
h_dim = get_h_dim(ds)
v_dim = get_v_dim(ds)
x0, y0, x1, y1 = (
(0, 1, 2, 3)
if (always_xy or is_xy_order(native_crs))
Expand All @@ -416,13 +421,13 @@ def _ensure_bbox_y_ascending(bbox: list, xy_order: bool = True):
bbox[y0], bbox[y1] = bbox[y1], bbox[y0]


def _get_h_dim(ds: xr.Dataset):
def get_h_dim(ds: xr.Dataset):
return [
d for d in list(map(str, ds.dims)) if d[:3].lower() in {'x', 'lon'}
][0]


def _get_v_dim(ds: xr.Dataset):
def get_v_dim(ds: xr.Dataset):
return [
d for d in list(map(str, ds.dims)) if d[:3].lower() in {'y', 'lat'}
][0]
Expand Down
122 changes: 122 additions & 0 deletions xcube/webapi/ows/coverages/scaling.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,122 @@
# The MIT License (MIT)
# Copyright (c) 2023 by the xcube team and contributors
#
# Permission is hereby granted, free of charge, to any person obtaining a
# copy of this software and associated documentation files (the "Software"),
# to deal in the Software without restriction, including without limitation
# the rights to use, copy, modify, merge, publish, distribute, sublicense,
# and/or sell copies of the Software, and to permit persons to whom the
# Software is furnished to do so, subject to the following conditions:
#
# The above copyright notice and this permission notice shall be included in
# all copies or substantial portions of the Software.
#
# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
# FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
# DEALINGS IN THE SOFTWARE.
from typing import Optional

import pyproj
import xarray as xr

from xcube.server.api import ApiError
from xcube.webapi.ows.coverages.controllers import get_h_dim, get_v_dim
from xcube.webapi.ows.coverages.request import CoveragesRequest


class CoverageScaling:

_scale: Optional[tuple[float, float]] = None
_final_size: Optional[tuple[int, int]] = None
_initial_size: tuple[int, int]
_crs: pyproj.CRS
_x_name: str
_y_name: str

def __init__(self, request: CoveragesRequest, crs: pyproj.CRS,
ds: xr.Dataset):
h_dim = get_h_dim(ds)
v_dim = get_v_dim(ds)
for d in h_dim, v_dim:
size = ds.dims[d]
if size == 0:
# Requirement 8C currently specifies a 204 rather than 404 here,
# but spec will soon be updated to allow 404 as an alternative.
# (J. Jacovella-St-Louis, pers. comm., 2023-11-27).
raise ApiError.NotFound(
f'Requested coverage contains no data: {d} has zero size.'
)
self._initial_size = ds.dims[h_dim], ds.dims[v_dim]

self._crs = crs
self._y_name = self.get_axis_from_crs(
{'lat', 'latitude', 'geodetic latitude', 'n', 'north', 'y'}
)
self._x_name = self.get_axis_from_crs(
{'longitude', 'geodetic longitude', 'lon', 'e', 'east', 'x'}
)

# The standard doesn't define behaviour when multiple scale
# parameters are given. We choose to handle one and ignore the
# others in such cases.
if request.scale_factor is not None:
self._scale = request.scale_factor, request.scale_factor
elif request.scale_axes is not None:
self._scale = self._get_xy_values(request.scale_axes)
elif request.scale_size is not None:
# The standard allows floats for "scale-size" but mandates:
# "The returned coverage SHALL contain exactly the specified number
# of sample values". We can't return a fractional number of sample
# values, so truncate to int here.
x, y = self._get_xy_values(request.scale_size)
self._final_size = int(x), int(y)
else:
# The standard doesn't mandate this as a default; in future, we
# might choose to downscale automatically if a large coverage
# is requested without an explicit scaling parameter.
self._scale = (1, 1)

@property
def scale(self) -> tuple[float, float]:
if self._scale is not None:
return self._scale
else:
x_initial, y_initial = self._initial_size
x_final, y_final = self._final_size
return x_initial / x_final, y_initial / y_final

@property
def size(self) -> tuple[float, float]:
if self._final_size is not None:
return self._final_size
else:
x_initial, y_initial = self._initial_size
x_scale, y_scale = self._scale
return x_initial / x_scale, y_initial / y_scale

def _get_xy_values(
self, axis_to_value: dict[str, float]
) -> tuple[float, float]:
x, y = None, None
for axis in axis_to_value:
if axis.lower()[:3] in ['x', 'e', 'eas', 'lon', self._x_name]:
x = axis_to_value[axis]
if axis.lower()[:3] in ['y', 'n', 'nor', 'lat', self._y_name]:
y = axis_to_value[axis]
return x, y

def get_axis_from_crs(self, valid_identifiers: set[str]):
for axis in self._crs.axis_info:
if not hasattr(axis, 'abbrev'):
continue

Check warning on line 115 in xcube/webapi/ows/coverages/scaling.py

View check run for this annotation

Codecov / codecov/patch

xcube/webapi/ows/coverages/scaling.py#L115

Added line #L115 was not covered by tests
identifiers = set(map(
lambda attr: getattr(axis, attr, '').lower(),
['name', 'abbrev']
))
if identifiers & valid_identifiers:
return axis.abbrev
return None

Check warning on line 122 in xcube/webapi/ows/coverages/scaling.py

View check run for this annotation

Codecov / codecov/patch

xcube/webapi/ows/coverages/scaling.py#L122

Added line #L122 was not covered by tests

0 comments on commit d47f59e

Please sign in to comment.