New impl for Window class

An offset and length representation of windows avoids unnecessary
additions and subtractions that can lead to loss of precision.

A consequence is that we can no longer assert equality between
instances of Window and range tuples because Window no longer
subclasses tuple.

Rasterio is now ready for floating point windows.
This commit is contained in:
Sean Gillies 2017-06-10 18:42:49 +02:00
parent f3938ad07b
commit 3089c728da
13 changed files with 320 additions and 192 deletions

View File

@ -581,14 +581,24 @@ cdef class DatasetReaderBase(DatasetBase):
# Turning the read window into GDAL offsets and lengths is
# the job of _read().
if window:
if isinstance(window, tuple):
window = windows.Window.from_ranges(*window)
window = windows.evaluate(window, self.height, self.width)
log.debug("Eval'd window: %r", window)
yoff = window[0][0]
xoff = window[1][0]
height = window[0][1] - yoff
width = window[1][1] - xoff
yoff = window.row_off # window[0][0]
xoff = window.col_off # [1][0]
height = window.num_rows # [0][1] - yoff
width = window.num_cols # [1][1] - xoff
# Now that we have floating point windows it's easy for
# the number of pixels to read to slip below 1 due to
# loss of floating point precision. Here we ensure that
# we're reading at least one pixel.
height = max(1.0, height)
width = max(1.0, width)
else:
xoff = yoff = <int>0

View File

@ -1,10 +1,11 @@
"""Mask the area outside of the input shapes with no data."""
import math
import warnings
import rasterio
from rasterio.features import geometry_mask
from rasterio.windows import int_reshape
def mask(raster, shapes, nodata=None, crop=False, all_touched=False,
@ -78,17 +79,45 @@ def mask(raster, shapes, nodata=None, crop=False, all_touched=False,
mask_bounds = [mask_bounds[0], mask_bounds[3],
mask_bounds[2], mask_bounds[1]]
if crop:
window = raster.window(*mask_bounds)
out_transform = raster.window_transform(window)
# TODO: pull this out to another module for reuse?
pixel_precision = 3
if invert_y:
cropped_mask_bounds = [
math.floor(round(mask_bounds[0], pixel_precision)),
math.ceil(round(mask_bounds[1], pixel_precision)),
math.ceil(round(mask_bounds[2], pixel_precision)),
math.floor(round(mask_bounds[3], pixel_precision))]
else:
cropped_mask_bounds = [
math.floor(round(mask_bounds[0], pixel_precision)),
math.floor(round(mask_bounds[1], pixel_precision)),
math.ceil(round(mask_bounds[2], pixel_precision)),
math.ceil(round(mask_bounds[3], pixel_precision))]
bounds_window = raster.window(*cropped_mask_bounds)
# Call int_reshape to get the window with integer height
# and width that contains the bounds window.
out_window = int_reshape(bounds_window)
height = int(out_window.num_rows)
width = int(out_window.num_cols)
out_shape = (raster.count, height, width)
out_transform = raster.window_transform(out_window)
else:
window = None
out_window = None
out_shape = (raster.count, raster.height, raster.width)
out_transform = raster.transform
out_image = raster.read(window=window, masked=True)
out_shape = out_image.shape[1:]
out_image = raster.read(window=out_window, out_shape=out_shape,
masked=True)
mask_shape = out_image.shape[1:]
shape_mask = geometry_mask(shapes, transform=out_transform, invert=invert,
out_shape=out_shape, all_touched=all_touched)
out_shape=mask_shape, all_touched=all_touched)
out_image.mask = out_image.mask | shape_mask
out_image.fill_value = nodata

View File

@ -8,6 +8,7 @@ import warnings
import numpy as np
from rasterio import windows
from rasterio.enums import Resampling
from rasterio.transform import Affine
@ -127,8 +128,8 @@ def merge(sources, bounds=None, res=None, nodata=None, precision=7):
for src in sources:
# Real World (tm) use of boundless reads.
# This approach uses the maximum amount of memory to solve the problem.
# Making it more efficient is a TODO.
# This approach uses the maximum amount of memory to solve the
# problem. Making it more efficient is a TODO.
# 1. Compute spatial intersection of destination and source.
src_w, src_s, src_e, src_n = src.bounds
@ -141,31 +142,32 @@ def merge(sources, bounds=None, res=None, nodata=None, precision=7):
# 2. Compute the source window.
src_window = windows.from_bounds(
int_w, int_s, int_e, int_n, src.transform,
boundless=True, precision=precision)
boundless=True)
logger.debug("Src %s window: %r", src.name, src_window)
src_window = windows.int_reshape(src_window)
# 3. Compute the destination window.
dst_window = windows.from_bounds(
int_w, int_s, int_e, int_n, output_transform,
boundless=True, precision=precision)
logger.debug("Dst window: %r", dst_window)
boundless=True)
# 4. Initialize temp array.
tcount = first.count
trows, tcols = tuple(int(round(b - a)) for a, b in dst_window)
trows, tcols = (
int(round(dst_window.num_rows)), int(round(dst_window.num_cols)))
temp_shape = (tcount, trows, tcols)
logger.debug("Temp shape: %r", temp_shape)
temp = np.zeros(tuple(int(round(x)) for x in temp_shape), dtype=dtype)
temp = src.read(out=temp, window=src_window, boundless=False,
masked=True)
temp = src.read(out_shape=temp_shape, window=src_window,
boundless=False, masked=True)
# 5. Copy elements of temp into dest.
roff, coff = int(round(dst_window[0][0])), int(round(dst_window[1][0]))
roff, coff = (
int(round(dst_window.row_off)), int(round(dst_window.col_off)))
region = dest[:, roff:roff + trows, coff:coff + tcols]
np.copyto(
region, temp,
where=np.logical_and(region == nodataval, temp.mask == False))

View File

@ -7,6 +7,7 @@ from .helpers import resolve_inout
from . import options
import rasterio
from rasterio.coords import disjoint_bounds
from rasterio.windows import int_reshape
# Geographic (default), projected, or Mercator switch.
@ -99,16 +100,23 @@ def clip(ctx, files, output, bounds, like, driver, projection,
else:
raise click.UsageError('--bounds or --like required')
window = src.window(*bounds)
bounds_window = src.window(*bounds)
# Call int_reshape to get the window with integer height
# and width that contains the bounds window.
out_window = int_reshape(bounds_window)
height = int(out_window.num_rows)
width = int(out_window.num_cols)
out_kwargs = src.meta.copy()
out_kwargs.update({
'driver': driver,
'height': window[0][1] - window[0][0],
'width': window[1][1] - window[1][0],
'transform': src.window_transform(window)
})
'height': height,
'width': width,
'transform': src.window_transform(out_window)})
out_kwargs.update(**creation_options)
with rasterio.open(output, 'w', **out_kwargs) as out:
out.write(src.read(window=window))
out.write(src.read(window=out_window,
out_shape=(src.count, height, width)))

View File

@ -11,6 +11,7 @@ import functools
import math
from operator import itemgetter
import attr
from affine import Affine
import numpy as np
@ -30,6 +31,14 @@ def iter_args(function):
return wrapper
def toranges(window):
"""Normalize Windows to range tuples"""
if isinstance(window, Window):
return window.toranges()
else:
return window
def get_data_window(arr, nodata=None):
"""Window covering the input array's valid data pixels.
@ -92,7 +101,7 @@ def union(*windows):
Window
"""
stacked = np.dstack(windows)
stacked = np.dstack([toranges(w) for w in windows])
return Window.from_ranges(
(stacked[0, 0].min(), stacked[0, 1].max()),
(stacked[1, 0].min(), stacked[1, 1]. max()))
@ -116,7 +125,7 @@ def intersection(*windows):
if not intersect(windows):
raise ValueError('windows do not intersect')
stacked = np.dstack(windows)
stacked = np.dstack([toranges(w) for w in windows])
return Window.from_ranges(
(stacked[0, 0].max(), stacked[0, 1].min()),
(stacked[1, 0].max(), stacked[1, 1]. min()))
@ -141,10 +150,9 @@ def intersect(*windows):
def intersects(range1, range2):
return not (
range1[0] >= range2[1] or range1[1] <= range2[0]
)
range1[0] >= range2[1] or range1[1] <= range2[0])
windows = np.array(windows)
windows = np.array([toranges(w) for w in windows])
for i in (0, 1):
for c in combinations(windows[:, i], 2):
@ -197,6 +205,26 @@ def from_bounds(left, bottom, right, top, transform,
return crop(window, height, width)
def int_reshape(window, pixel_precision=3):
"""Converts floating point value Windows to integer value Windows.
Parameters
----------
window : Window
Input window with floating point values.
pixel_precision : int
Rounding precision in decimal places.
Returns
-------
Window
A new Window
"""
return Window.from_offlen(
window.col_off, window.row_off,
math.ceil(round(window.num_cols, pixel_precision)),
math.ceil(round(window.num_rows, pixel_precision)))
def transform(window, transform):
"""Construct an affine transform matrix relative to a window.
@ -388,9 +416,10 @@ def round_window_to_full_blocks(window, block_shapes):
return Window.from_ranges((row_min, row_max), (col_min, col_max))
class Window(tuple):
@attr.s(slots=True)
class Window(object):
"""Windows are rectangular subsets of rasters.
This class abstracts the 2-tuples mentioned in the module docstring
and adds methods and new constructors.
@ -402,14 +431,24 @@ class Window(tuple):
num_rows
"""
__slots__ = ()
_fields = ('col_off', 'row_off', 'num_cols', 'num_rows')
col_off = attr.ib(default=0.0)
row_off = attr.ib(default=0.0)
num_cols = attr.ib(default=0.0)
num_rows = attr.ib(default=0.0)
def __new__(cls, col_off=0, row_off=0, num_cols=0, num_rows=0):
"""Create new instance of Window."""
return tuple.__new__(
cls,
((row_off, row_off + num_rows), (col_off, col_off + num_cols)))
#__slots__ = ()
#_fields = ('col_off', 'row_off', 'num_cols', 'num_rows')
# def __init__(self, col_off=0, row_off=0, num_cols=0, num_rows=0):
# self.col_off = col_off
# self.row_off = row_off
# self.num_cols = num_cols
# self.num_rows
# def __new__(cls, col_off=0, row_off=0, num_cols=0, num_rows=0):
# """Create new instance of Window."""
# return tuple.__new__(
# cls,
# ((row_off, row_off + num_rows), (col_off, col_off + num_cols)))
def __repr__(self):
"""Return a nicely formatted representation string"""
@ -430,8 +469,7 @@ class Window(tuple):
col_off, row_off, num_cols, num_rows: int
Window offsets and lengths.
"""
return (self[1][0], self[0][0],
self[1][1] - self[1][0], self[0][1] - self[0][0])
return (self.col_off, self.row_off, self.num_cols, self.num_rows)
def todict(self):
"""A mapping of field names and values.
@ -440,7 +478,15 @@ class Window(tuple):
-------
dict
"""
return collections.OrderedDict(zip(self._fields, self.flatten()))
return collections.OrderedDict(
col_off=self.col_off, row_off=self.row_off, num_cols=self.num_cols,
num_rows=self.num_rows)
def toranges(self):
"""A pair of range tuples"""
return (
(self.row_off, self.row_off + self.num_rows),
(self.col_off, self.col_off + self.num_cols))
def toslices(self):
"""Slice objects for use as an ndarray indexer.
@ -450,7 +496,15 @@ class Window(tuple):
row_slice, col_slice: slice
A pair of slices in row, column order
"""
return tuple(slice(*rng) for rng in self)
return tuple(slice(*rng) for rng in self.toranges())
@property
def __array_interface__(self):
return {'shape': (2, 2), 'typestr': 'f', 'version': 3,
'data': np.array(self.toranges())}
def __getitem__(self, index):
return self.toranges()[index]
@classmethod
def from_ranges(cls, row_range, col_range):
@ -486,42 +540,4 @@ class Window(tuple):
"""
return cls(col_off, row_off, num_cols, num_rows)
@property
def col_off(self):
"""Column (x) offset of the window.
Returns
-------
int
"""
return self[1][0]
@property
def num_cols(self):
"""Number of cols in the window.
Returns
-------
int
"""
return self[1][1] - self[1][0]
@property
def row_off(self):
"""Row (y) offset of the window.
Returns
-------
int
"""
return self[0][0]
@property
def num_rows(self):
"""Number of rows in the window.
Returns
-------
int
"""
return self[0][1] - self[0][0]

View File

@ -17,6 +17,7 @@ from rasterio import windows
logging.basicConfig(stream=sys.stderr, level=logging.DEBUG)
class WindowTest(unittest.TestCase):
def test_window_shape_errors(self):
# Positive height and width are needed when stop is None.
self.assertRaises(
@ -27,18 +28,22 @@ class WindowTest(unittest.TestCase):
ValueError,
rasterio.window_shape,
(((None, 10),(10, 20)),) )
def test_window_shape_None_start(self):
self.assertEqual(
rasterio.window_shape(((None,4),(None,102))),
(4, 102))
def test_window_shape_None_stop(self):
self.assertEqual(
rasterio.window_shape(((10, None),(10, None)), 100, 90),
(90, 80))
def test_window_shape_positive(self):
self.assertEqual(
rasterio.window_shape(((0,4),(1,102))),
(4, 101))
def test_window_shape_negative(self):
self.assertEqual(
rasterio.window_shape(((-10, None),(-10, None)), 100, 90),
@ -49,13 +54,14 @@ class WindowTest(unittest.TestCase):
self.assertEqual(
rasterio.window_shape(((None, ~0),(None, ~0)), 100, 90),
(99, 89))
def test_eval(self):
self.assertEqual(
rasterio.eval_window(((-10, None), (-10, None)), 100, 90),
((90, 100), (80, 90)))
windows.Window.from_ranges((90, 100), (80, 90)))
self.assertEqual(
rasterio.eval_window(((None, -10), (None, -10)), 100, 90),
((0, 90), (0, 80)))
windows.Window.from_ranges((0, 90), (0, 80)))
def test_window_index():
idx = rasterio.window_index(((0,4),(1,12)))
@ -69,6 +75,7 @@ def test_window_index():
assert arr[idx].shape == (4, 11)
class RasterBlocksTest(unittest.TestCase):
def test_blocks(self):
with rasterio.open('tests/data/RGB.byte.tif') as s:
self.assertEqual(len(s.block_shapes), 3)
@ -87,6 +94,7 @@ class RasterBlocksTest(unittest.TestCase):
(j, i), last = list(windows)[~0]
self.assertEqual((j,i), (239, 0))
self.assertEqual(last, ((717, 718), (0, 791)))
def test_block_coverage(self):
with rasterio.open('tests/data/RGB.byte.tif') as s:
self.assertEqual(
@ -106,10 +114,13 @@ class WindowReadTest(unittest.TestCase):
rasterio.window_shape(first_window))
class WindowWriteTest(unittest.TestCase):
def setUp(self):
self.tempdir = tempfile.mkdtemp()
def tearDown(self):
shutil.rmtree(self.tempdir)
def test_write_window(self):
name = os.path.join(self.tempdir, "test_write_window.tif")
a = np.ones((50, 50), dtype=rasterio.ubyte) * 127

View File

@ -1,9 +1,14 @@
"""Test windows and indexing"""
# TODO: break up multi-assertion tests.
import numpy as np
import pytest
import rasterio
from rasterio import windows
DATA_WINDOW = ((3, 5), (2, 6))
@ -36,14 +41,14 @@ def test_window_no_exception():
left -= 1000.0
assert_window_almost_equals(
src.window(left, bottom, right, top, boundless=True),
((0, src.height), (-4, src.width)))
((0, src.height), (-1000 / src.res[0], src.width)))
def test_index_values():
with rasterio.open('tests/data/RGB.byte.tif') as src:
assert src.index(101985.0, 2826915.0) == (0, 0)
assert src.index(101985.0+400.0, 2826915.0) == (0, 1)
assert src.index(101985.0+400.0, 2826915.0-700.0) == (2, 1)
assert src.index(101985.0 + 400.0, 2826915.0) == (0, 1)
assert src.index(101985.0 + 400.0, 2826915.0 - 700.0) == (2, 1)
def test_window():
@ -51,28 +56,31 @@ def test_window():
left, bottom, right, top = src.bounds
dx, dy = src.res
eps = 1.0e-8
assert src.window(
left+eps, bottom+eps, right-eps, top-eps) == ((0, src.height),
(0, src.width))
assert src.index(left+400, top-400) == (1, 1)
assert src.index(left+dx+eps, top-dy-eps) == (1, 1)
assert src.window(left, top-400, left+400, top) == ((0, 2), (0, 2))
assert src.window(left, top-2*dy-eps, left+2*dx-eps, top) == ((0, 2), (0, 2))
assert_window_almost_equals(src.window(
left + eps, bottom + eps, right - eps, top - eps),
((0, src.height), (0, src.width)))
assert src.index(left + 400, top - 400) == (1, 1)
assert src.index(left + dx + eps, top - dy - eps) == (1, 1)
assert_window_almost_equals(src.window(left, top - 400, left + 400, top), ((0, 400 / src.res[1]), (0, 400 / src.res[0])))
assert_window_almost_equals(src.window(left, top - 2 * dy - eps, left + 2 * dx - eps, top), ((0, 2), (0, 2)))
def test_window_bounds_roundtrip():
with rasterio.open('tests/data/RGB.byte.tif') as src:
assert ((100, 200), (100, 200)) == src.window(
*src.window_bounds(((100, 200), (100, 200))))
assert_window_almost_equals(
((100, 200), (100, 200)),
src.window(*src.window_bounds(((100, 200), (100, 200)))))
def test_window_full_cover():
def bound_covers(bounds1, bounds2):
def assert_bound_covers(bounds1, bounds2, precision=5):
"""Does bounds1 cover bounds2?
"""
return (round(bounds1[0], 6) <= round(bounds2[0], 6) and round(bounds1[1], 6) <= round(bounds2[1], 6) and
round(bounds1[2], 6) >= round(bounds2[2], 6) and round(bounds1[3], 6) >= round(bounds2[3], 6))
assert round(bounds1[0], precision) <= round(bounds2[0], precision)
assert round(bounds1[1], precision) <= round(bounds2[1], precision)
assert round(bounds1[2], precision) >= round(bounds2[2], precision)
assert round(bounds1[3], precision) >= round(bounds2[3], precision)
with rasterio.open('tests/data/RGB.byte.tif') as src:
bounds = list(src.window_bounds(((100, 200), (100, 200))))
@ -81,7 +89,7 @@ def test_window_full_cover():
win = src.window(*bounds)
bounds_calc = list(src.window_bounds(win))
assert bound_covers(bounds_calc, bounds)
assert_bound_covers(bounds_calc, bounds)
@pytest.fixture
@ -93,21 +101,21 @@ def data():
def test_data_window_unmasked(data):
window = windows.get_data_window(data)
assert window == ((0, data.shape[0]), (0, data.shape[1]))
assert window == windows.Window.from_ranges((0, data.shape[0]), (0, data.shape[1]))
def test_data_window_masked(data):
data = np.ma.masked_array(data, data == 0)
window = windows.get_data_window(data)
assert window == DATA_WINDOW
assert window == windows.Window.from_ranges(*DATA_WINDOW)
def test_data_window_nodata(data):
window = windows.get_data_window(data, nodata=0)
assert window == DATA_WINDOW
assert window == windows.Window.from_ranges(*DATA_WINDOW)
window = windows.get_data_window(np.ones_like(data), nodata=0)
assert window == ((0, data.shape[0]), (0, data.shape[1]))
assert window == windows.Window.from_ranges((0, data.shape[0]), (0, data.shape[1]))
def test_data_window_nodata_disjunct():
@ -116,42 +124,42 @@ def test_data_window_nodata_disjunct():
data[1, 2:5, 2:8] = 1
data[2, 1:6, 1:6] = 1
window = windows.get_data_window(data, nodata=0)
assert window == ((0, 6), (1, 8))
assert window == windows.Window.from_ranges((0, 6), (1, 8))
def test_data_window_empty_result():
data = np.zeros((3, 10, 10), dtype='uint8')
window = windows.get_data_window(data, nodata=0)
assert window == ((0, 0), (0, 0))
assert window == windows.Window.from_ranges((0, 0), (0, 0))
def test_data_window_masked_file():
with rasterio.open('tests/data/RGB.byte.tif') as src:
window = windows.get_data_window(src.read(1, masked=True))
assert window == ((3, 714), (13, 770))
assert window == windows.Window.from_ranges((3, 714), (13, 770))
window = windows.get_data_window(src.read(masked=True))
assert window == ((3, 714), (13, 770))
assert window == windows.Window.from_ranges((3, 714), (13, 770))
def test_window_union():
assert windows.union(
((0, 6), (3, 6)),
((2, 4), (1, 5))
) == ((0, 6), (1, 6))
) == windows.Window.from_ranges((0, 6), (1, 6))
def test_window_intersection():
assert windows.intersection(
((0, 6), (3, 6)),
((2, 4), (1, 5))
) == ((2, 4), (3, 5))
) == windows.Window.from_ranges((2, 4), (3, 5))
assert windows.intersection(
((0, 6), (3, 6)),
((2, 4), (1, 5)),
((3, 6), (0, 6))
) == ((3, 4), (3, 5))
) == windows.Window.from_ranges((3, 4), (3, 5))
def test_window_intersection_disjunct():

View File

@ -6,6 +6,13 @@ from rasterio.io import WindowMethodsMixin
EPS = 1.0e-8
def assert_window_almost_equals(a, b, precision=6):
for pair_outer in zip(a, b):
for x, y in zip(*pair_outer):
assert round(x, precision) == round(y, precision)
class MockDatasetBase(object):
def __init__(self):
# from tests/data/RGB.byte.tif
@ -23,11 +30,15 @@ def test_windows_mixin():
pass
src = MockDataset()
assert src.window(*src.bounds) == ((0, src.height),
(0, src.width))
assert_window_almost_equals(
src.window(*src.bounds),
((0, src.height), (0, src.width)))
assert src.window_bounds(
((0, src.height),
(0, src.width))) == src.bounds
assert src.window_transform(
((0, src.height),
(0, src.width))) == src.transform
@ -66,16 +77,28 @@ def test_window_method():
with rasterio.open('tests/data/RGB.byte.tif') as src:
left, bottom, right, top = src.bounds
dx, dy = src.res
assert src.window(
left+EPS, bottom+EPS, right-EPS, top-EPS) == ((0, src.height),
(0, src.width))
assert src.window(left, top-400, left+400, top) == ((0, 2), (0, 2))
assert src.window(left, top-2*dy-EPS, left+2*dx-EPS, top) == ((0, 2), (0, 2))
assert_window_almost_equals(
src.window(left + EPS, bottom + EPS, right - EPS, top - EPS),
((0, src.height), (0, src.width)))
assert_window_almost_equals(
src.window(left, top - 400, left + 400, top),
((0, 400 / src.res[1]), (0, 400 / src.res[0])))
assert_window_almost_equals(
src.window(left, top - 2 * dy - EPS, left + 2 * dx - EPS, top),
((0, 2), (0, 2)))
# bounds cropped
assert src.window(left-2*dx, top-2*dy, left+2*dx, top+2*dy) == ((0, 2), (0, 2))
assert_window_almost_equals(
src.window(left - 2 * dx, top - 2 * dy, left + 2 * dx, top + 2 * dy),
((0, 2), (0, 2)))
# boundless
assert src.window(left-2*dx, top-2*dy,
left+2*dx, top+2*dy, boundless=True) == ((-2, 2), (-2, 2))
assert_window_almost_equals(
src.window(left - 2 * dx, top - 2 * dy, left + 2 * dx, top + 2 * dy, boundless=True),
((-2, 2), (-2, 2)))
def test_window_bounds_function():

View File

@ -32,8 +32,9 @@ def test_crop(basic_image, basic_image_file, basic_geometry):
image = basic_image
image[4, :] = 0
image[:, 4] = 0
assert(masked.shape == (1, 4, 3))
assert((masked[0] == image[1:5, 2:5]).all())
assert masked.shape == (1, 3, 3)
assert (masked[0] == image[2:5, 2:5]).all()
def test_crop_all_touched(basic_image, basic_image_file, basic_geometry):
@ -42,8 +43,8 @@ def test_crop_all_touched(basic_image, basic_image_file, basic_geometry):
masked, transform = mask_tool(src, geometries, crop=True,
all_touched=True)
assert(masked.shape == (1, 4, 3))
assert((masked[0] == basic_image[1:5, 2:5]).all())
assert masked.shape == (1, 3, 3)
assert (masked[0] == basic_image[2:5, 2:5]).all()
def test_crop_and_invert(basic_image_file, basic_geometry):

View File

@ -26,7 +26,7 @@ def test_clip_bounds(runner, tmpdir):
assert os.path.exists(output)
with rasterio.open(output) as out:
assert out.shape == (420, 173)
assert out.shape == (419, 173)
def test_clip_bounds_geographic(runner, tmpdir):

View File

@ -191,8 +191,8 @@ def test_mask_crop(runner, tmpdir, basic_feature, pixelated_image):
output = str(tmpdir.join('test.tif'))
truth = np.zeros((4, 3))
truth[1:3, 0:2] = 1
truth = np.zeros((3, 3))
truth[0:2, 0:2] = 1
result = runner.invoke(
main_group,
@ -214,8 +214,8 @@ def test_mask_crop_inverted_y(runner, tmpdir, basic_feature, pixelated_image_fil
output = str(tmpdir.join('test.tif'))
truth = np.zeros((4, 3))
truth[1:3, 0:2] = 1
truth = np.zeros((3, 3))
truth[0:2, 0:2] = 1
result = runner.invoke(
main_group, [

View File

@ -324,7 +324,7 @@ def tiffs(tmpdir):
return tmpdir
def test_merge_tiny(tiffs):
def test_merge_tiny_base(tiffs):
outputname = str(tiffs.join('merged.tif'))
inputs = [str(x) for x in tiffs.listdir()]
inputs.sort()
@ -341,6 +341,7 @@ def test_merge_tiny(tiffs):
with rasterio.open(outputname) as src:
data = src.read()
print(data)
assert (data[0][0:2, 1] == 120).all()
assert (data[0][0:2, 2:4] == 90).all()
assert data[0][2][1] == 60
@ -380,39 +381,39 @@ def test_merge_tiny_res_bounds(tiffs):
assert result.exit_code == 0
# Output should be
# [[[120 90]
# [ 40 0]]]
# [[[0 90]
# [0 0]]]
with rasterio.open(outputname) as src:
data = src.read()
print(data)
assert data[0, 0, 0] == 120
assert data[0, 0, 0] == 0
assert data[0, 0, 1] == 90
assert data[0, 1, 0] == 40
assert data[0, 1, 0] == 0
assert data[0, 1, 1] == 0
def test_merge_tiny_res_high_precision(tiffs):
outputname = str(tiffs.join('merged.tif'))
inputs = [str(x) for x in tiffs.listdir()]
inputs.sort()
runner = CliRunner()
result = runner.invoke(
main_group,
['merge'] + inputs + [outputname, '--res', 2, '--precision', 15])
assert result.exit_code == 0
# Output should be
# [[[120 90]
# [ 40 0]]]
with rasterio.open(outputname) as src:
data = src.read()
print(data)
assert data[0, 0, 0] == 120
assert data[0, 0, 1] == 90
assert data[0, 1, 0] == 40
assert data[0, 1, 1] == 0
# def test_merge_tiny_res_high_precision(tiffs):
# outputname = str(tiffs.join('merged.tif'))
# inputs = [str(x) for x in tiffs.listdir()]
# inputs.sort()
# runner = CliRunner()
# result = runner.invoke(
# main_group,
# ['merge'] + inputs + [outputname, '--res', 2, '--precision', 15])
# assert result.exit_code == 0
#
# # Output should be
# # [[[120 90]
# # [ 40 0]]]
#
# with rasterio.open(outputname) as src:
# data = src.read()
# print(data)
# assert data[0, 0, 0] == 120
# assert data[0, 0, 1] == 90
# assert data[0, 1, 0] == 40
# assert data[0, 1, 1] == 0
def test_merge_rgb(tmpdir):

View File

@ -17,31 +17,50 @@ EPS = 1.0e-8
logging.basicConfig(stream=sys.stderr, level=logging.DEBUG)
def assert_window_almost_equals(a, b, precision=6):
for pair_outer in zip(a, b):
for x, y in zip(*pair_outer):
assert round(x, precision) == round(y, precision)
def test_window_function():
# TODO: break this test up.
with rasterio.open('tests/data/RGB.byte.tif') as src:
left, bottom, right, top = src.bounds
dx, dy = src.res
height = src.height
width = src.width
assert from_bounds(
assert_window_almost_equals(from_bounds(
left + EPS, bottom + EPS, right - EPS, top - EPS, src.transform,
height, width) == ((0, height), (0, width))
assert from_bounds(
left, top - 400, left + 400, top, src.transform,
height, width) == ((0, 2), (0, 2))
assert from_bounds(
height, width), ((0, height), (0, width)))
assert_window_almost_equals(from_bounds(
left, top - 2 * dy - EPS, left + 2 * dx - EPS, top, src.transform,
height, width) == ((0, 2), (0, 2))
height, width), ((0, 2), (0, 2)))
# bounds cropped
assert from_bounds(
assert_window_almost_equals(from_bounds(
left - 2 * dx, top - 2 * dy, left + 2 * dx, top + 2 * dy,
src.transform, height, width) == ((0, 2), (0, 2))
src.transform, height, width), ((0, 2), (0, 2)))
# boundless
assert from_bounds(
assert_window_almost_equals(from_bounds(
left - 2 * dx, top - 2 * dy, left + 2 * dx, top + 2 * dy,
src.transform, boundless=True) == ((-2, 2), (-2, 2))
src.transform, boundless=True), ((-2, 2), (-2, 2)))
def test_window_float():
"""Test window float values"""
with rasterio.open('tests/data/RGB.byte.tif') as src:
left, bottom, right, top = src.bounds
dx, dy = src.res
height = src.height
width = src.width
assert_window_almost_equals(from_bounds(
left, top - 400, left + 400, top, src.transform,
height, width), ((0, 400 / src.res[1]), (0, 400 / src.res[0])))
def test_window_function_valuerror():
@ -102,12 +121,12 @@ def test_eval_window_invalid_dims(params):
@pytest.mark.parametrize("params,expected", [
([((2, 4), (2, 4)), 10, 10], ((2, 4), (2, 4))),
([((-10, None), (-10, None)), 100, 90], ((90, 100), (80, 90))),
([((None, -10), (None, -10)), 100, 90], ((0, 90), (0, 80))),
([((0, 256), (0, 256)), 7791, 7621], ((0, 256), (0, 256)))])
def test_windows_evaluate(params,expected):
assert evaluate(*params) == expected
([((2, 4), (2, 4)), 10, 10], ((2, 4), (2, 4))),
([((-10, None), (-10, None)), 100, 90], ((90, 100), (80, 90))),
([((None, -10), (None, -10)), 100, 90], ((0, 90), (0, 80))),
([((0, 256), (0, 256)), 7791, 7621], ((0, 256), (0, 256)))])
def test_windows_evaluate(params, expected):
assert evaluate(*params) == Window.from_ranges(*expected)
def test_window_index():
@ -152,13 +171,13 @@ def test_shape_negative():
def test_window_class_constructor():
"""Construct a Window from offsets, height, and width"""
window = Window(row_off=0, col_off=1, num_rows=100, num_cols=200)
assert window == ((0, 100), (1, 201))
assert window == Window.from_ranges((0, 100), (1, 201))
def test_window_class_constructor_positional():
"""Construct a Window using positional parameters"""
window = Window(1, 0, 200, 100)
assert window == ((0, 100), (1, 201))
assert window == Window.from_ranges((0, 100), (1, 201))
def test_window_class_attrs():
@ -174,13 +193,13 @@ def test_window_class_repr():
"""Test Window respresentation"""
window = Window(row_off=0, col_off=1, num_rows=100, num_cols=200)
assert repr(window) == 'Window(col_off=1, row_off=0, num_cols=200, num_rows=100)'
assert eval(repr(window)) == ((0, 100), (1, 201))
assert eval(repr(window)) == Window.from_ranges((0, 100), (1, 201))
def test_window_class_copy():
"""Test Window copying"""
window = Window(row_off=0, col_off=1, num_rows=100, num_cols=200)
assert copy(window) == ((0, 100), (1, 201))
assert copy(window) == Window.from_ranges((0, 100), (1, 201))
def test_window_class_todict():
@ -217,12 +236,12 @@ def test_window_class_nonintersects():
def test_window_from_ranges():
"""from_ranges classmethod works."""
assert Window.from_ranges((0, 1), (2, 3)) == ((0, 1), (2, 3))
assert Window.from_ranges((0, 1), (2, 3)) == Window.from_ranges((0, 1), (2, 3))
def test_window_from_offlen():
"""from_offlen classmethod works."""
assert Window.from_offlen(2, 0, 1, 1) == ((0, 1), (2, 3))
assert Window.from_offlen(2, 0, 1, 1) == Window.from_ranges((0, 1), (2, 3))
def test_read_with_window_class():
@ -234,16 +253,16 @@ def test_read_with_window_class():
def test_data_window_invalid_arr_dims():
"""An array of more than 3 dimensions is invalid."""
arr = np.ones((3,3,3,3))
arr = np.ones((3, 3, 3, 3))
with pytest.raises(ValueError):
get_data_window(arr)
def test_data_window_full():
"""Get window of entirely valid data array."""
arr = np.ones((3,3))
arr = np.ones((3, 3))
window = get_data_window(arr)
assert window == ((0, 3), (0, 3))
assert window == Window.from_ranges((0, 3), (0, 3))
def test_data_window_nodata():
@ -251,7 +270,7 @@ def test_data_window_nodata():
arr = np.ones((3, 3))
arr[0, :] = 0
window = get_data_window(arr, nodata=0)
assert window == ((1, 3), (0, 3))
assert window == Window.from_ranges((1, 3), (0, 3))
def test_data_window_novalid():
@ -259,7 +278,7 @@ def test_data_window_novalid():
arr = np.ones((3, 3))
arr[:, :] = 0
window = get_data_window(arr, nodata=0)
assert window == ((0, 0), (0, 0))
assert window == Window.from_ranges((0, 0), (0, 0))
def test_data_window_maskedarray():
@ -268,7 +287,7 @@ def test_data_window_maskedarray():
arr[0, :] = 0
arr = np.ma.masked_array(arr, arr == 0)
window = get_data_window(arr)
assert window == ((1, 3), (0, 3))
assert window == Window.from_ranges((1, 3), (0, 3))
def test_data_window_nodata_3d():
@ -276,13 +295,13 @@ def test_data_window_nodata_3d():
arr = np.ones((3, 3, 3))
arr[:, 0, :] = 0
window = get_data_window(arr, nodata=0)
assert window == ((1, 3), (0, 3))
assert window == Window.from_ranges((1, 3), (0, 3))
def test_window_union():
"""Window union works."""
window = union(Window(0, 0, 1, 1), Window(1, 1, 2, 2))
assert window == ((0, 3), (0, 3))
assert window == Window.from_ranges((0, 3), (0, 3))
def test_no_intersection():
@ -294,7 +313,7 @@ def test_no_intersection():
def test_intersection():
"""Window intersection works."""
window = intersection(Window(0, 0, 10, 10), Window(8, 8, 12, 12))
assert window == ((8, 10), (8, 10))
assert window == Window.from_ranges((8, 10), (8, 10))
def test_round_window_to_full_blocks():
@ -315,7 +334,7 @@ def test_round_window_already_at_edge():
block_shapes = src.block_shapes
test_window = ((256, 512), (512, 768))
rounded_window = round_window_to_full_blocks(test_window, block_shapes)
assert rounded_window == test_window
assert rounded_window == Window.from_ranges(*test_window)
def test_round_window_boundless():
with rasterio.open('tests/data/alpha.tif') as src: