diff --git a/.coveragerc b/.coveragerc index d3a89f1d..d4a1f13a 100644 --- a/.coveragerc +++ b/.coveragerc @@ -2,9 +2,9 @@ plugins = Cython.Coverage source = rasterio omit = *.pxd +branch = True [report] -show_missing = True exclude_lines = pragma: no cover def __repr__ diff --git a/.gitignore b/.gitignore index 0b87ebdf..bd980805 100644 --- a/.gitignore +++ b/.gitignore @@ -42,6 +42,7 @@ htmlcov/ nosetests.xml coverage.xml *,cover +rasterio/_*.html # Translations *.mo diff --git a/rasterio/_features.pxd b/rasterio/_features.pxd index e4898ae0..1ca21ba1 100644 --- a/rasterio/_features.pxd +++ b/rasterio/_features.pxd @@ -12,7 +12,6 @@ cdef class GeomBuilder: cpdef _buildPolygon(self) cpdef _buildMultiPolygon(self) cdef build(self, void *geom) - cpdef build_wkb(self, object wkb) cdef class OGRGeomBuilder: diff --git a/rasterio/_features.pyx b/rasterio/_features.pyx index 02893030..89f9a3dc 100644 --- a/rasterio/_features.pyx +++ b/rasterio/_features.pyx @@ -1,6 +1,5 @@ import logging -import json import numpy as np cimport numpy as np from rasterio._io cimport InMemoryRaster @@ -62,10 +61,16 @@ def _shapes(image, mask, connectivity, transform): cdef InMemoryRaster mem_ds = None cdef InMemoryRaster mask_ds = None cdef bint is_float = np.dtype(image.dtype).kind == 'f' - cdef int fieldtp = 0 + cdef int fieldtp = 2 if is_float else 0 - if is_float: - fieldtp = 2 + valid_dtypes = ('int16', 'int32', 'uint8', 'uint16', 'float32') + + if np.dtype(image.dtype).name not in valid_dtypes: + raise ValueError('image dtype must be one of: %s' + % (', '.join(valid_dtypes))) + + if connectivity not in (4, 8): + raise ValueError("Connectivity Option must be 4 or 8") if dtypes.is_ndarray(image): mem_ds = InMemoryRaster(image, transform) @@ -76,17 +81,22 @@ def _shapes(image, mask, connectivity, transform): else: raise ValueError("Invalid source image") - if dtypes.is_ndarray(mask): - # A boolean mask must be converted to uint8 for GDAL - mask_ds = InMemoryRaster(mask.astype('uint8'), transform) - hmaskband = mask_ds.band - elif isinstance(mask, tuple): + if mask is not None: if mask.shape != image.shape: raise ValueError("Mask must have same shape as image") - mrdr = mask.ds - hmaskband = mrdr.band(mask.bidx) - else: - hmaskband = NULL + + if np.dtype(mask.dtype).name not in ('bool', 'uint8'): + raise ValueError("Mask must be dtype rasterio.bool_ or " + "rasterio.uint8") + + if dtypes.is_ndarray(mask): + # A boolean mask must be converted to uint8 for GDAL + mask_ds = InMemoryRaster(mask.astype('uint8'), transform) + hmaskband = mask_ds.band + + elif isinstance(mask, tuple): + mrdr = mask.ds + hmaskband = mrdr.band(mask.bidx) # Create an in-memory feature store. hfdriver = _ogr.OGRGetDriverByName("Memory") @@ -133,7 +143,7 @@ def _shapes(image, mask, connectivity, transform): _gdal.CSLDestroy(options) -def _sieve(image, size, output, mask, connectivity): +def _sieve(image, size, out, mask, connectivity): """ Replaces small polygons in `image` with the value of their largest neighbor. Polygons are found for each set of neighboring pixels of the @@ -147,7 +157,7 @@ def _sieve(image, size, output, mask, connectivity): rasterio.uint16, or rasterio.float32. size : int minimum polygon size (number of pixels) to retain. - output : numpy ndarray + out : numpy ndarray Array of same shape and data type as `image` in which to store results. mask : numpy ndarray or rasterio Band object Values of False or 0 will be excluded from feature generation. @@ -168,6 +178,29 @@ def _sieve(image, size, output, mask, connectivity): cdef _io.RasterUpdater udr cdef _io.RasterReader mask_reader + valid_dtypes = ('int16', 'int32', 'uint8', 'uint16') + + if np.dtype(image.dtype).name not in valid_dtypes: + valid_types_str = ', '.join(('rasterio.{0}'.format(t) for t + in valid_dtypes)) + raise ValueError('image dtype must be one of: %s' % valid_types_str) + + if size <= 0: + raise ValueError('size must be greater than 0') + elif type(size) == float: + raise ValueError('size must be an integer number of pixels') + elif size > (image.shape[0] * image.shape[1]): + raise ValueError('size must be smaller than size of image') + + if connectivity not in (4, 8): + raise ValueError('connectivity must be 4 or 8') + + if out.shape != image.shape: + raise ValueError('out raster shape must be same as image shape') + + if np.dtype(image.dtype).name != np.dtype(out.dtype).name: + raise ValueError('out raster must match dtype of image') + if dtypes.is_ndarray(image): in_mem_ds = InMemoryRaster(image) in_band = in_mem_ds.band @@ -177,27 +210,32 @@ def _sieve(image, size, output, mask, connectivity): else: raise ValueError("Invalid source image") - if dtypes.is_ndarray(output): - log.debug("Output array: %r", output) - out_mem_ds = InMemoryRaster(output) + if dtypes.is_ndarray(out): + log.debug("out array: %r", out) + out_mem_ds = InMemoryRaster(out) out_band = out_mem_ds.band - elif isinstance(output, tuple): - udr = output.ds - out_band = udr.band(output.bidx) + elif isinstance(out, tuple): + udr = out.ds + out_band = udr.band(out.bidx) else: - raise ValueError("Invalid output image") + raise ValueError("Invalid out image") - if dtypes.is_ndarray(mask): - # A boolean mask must be converted to uint8 for GDAL - mask_mem_ds = InMemoryRaster(mask.astype('uint8')) - mask_band = mask_mem_ds.band - elif isinstance(mask, tuple): + if mask is not None: if mask.shape != image.shape: raise ValueError("Mask must have same shape as image") - mask_reader = mask.ds - mask_band = mask_reader.band(mask.bidx) - else: - mask_band = NULL + + if np.dtype(mask.dtype) not in ('bool', 'uint8'): + raise ValueError("Mask must be dtype rasterio.bool_ or " + "rasterio.uint8") + + if dtypes.is_ndarray(mask): + # A boolean mask must be converted to uint8 for GDAL + mask_mem_ds = InMemoryRaster(mask.astype('uint8')) + mask_band = mask_mem_ds.band + + elif isinstance(mask, tuple): + mask_reader = mask.ds + mask_band = mask_reader.band(mask.bidx) _gdal.GDALSieveFilter( in_band, @@ -210,8 +248,8 @@ def _sieve(image, size, output, mask, connectivity): NULL ) - # Read from out_band into output - _io.io_auto(output, out_band, False) + # Read from out_band into out + _io.io_auto(out, out_band, False) if in_mem_ds is not None: in_mem_ds.close() @@ -238,7 +276,7 @@ def _rasterize(shapes, image, transform, all_touched): all_touched : boolean, optional If True, all pixels touched by geometries will be burned in. If false, only pixels whose center is within the polygon or - that are selected by brezenhams line algorithm will be burned + that are selected by Bresenham's line algorithm will be burned in. """ @@ -310,24 +348,21 @@ def _bounds(geometry): TODO: add to Fiona. """ - try: - if 'features' in geometry: - xmins = [] - ymins = [] - xmaxs = [] - ymaxs = [] - for feature in geometry['features']: - xmin, ymin, xmax, ymax = _bounds(feature['geometry']) - xmins.append(xmin) - ymins.append(ymin) - xmaxs.append(xmax) - ymaxs.append(ymax) - return min(xmins), min(ymins), max(xmaxs), max(ymaxs) - else: - xyz = tuple(zip(*list(_explode(geometry['coordinates'])))) - return min(xyz[0]), min(xyz[1]), max(xyz[0]), max(xyz[1]) - except (KeyError, TypeError): - return None + if 'features' in geometry: + xmins = [] + ymins = [] + xmaxs = [] + ymaxs = [] + for feature in geometry['features']: + xmin, ymin, xmax, ymax = _bounds(feature['geometry']) + xmins.append(xmin) + ymins.append(ymin) + xmaxs.append(xmax) + ymaxs.append(ymax) + return min(xmins), min(ymins), max(xmaxs), max(ymaxs) + else: + xyz = tuple(zip(*list(_explode(geometry['coordinates'])))) + return min(xyz[0]), min(xyz[1]), max(xyz[0]), max(xyz[1]) # Mapping of OGR integer geometry types to GeoJSON type names. @@ -359,16 +394,6 @@ GEOJSON2OGR_GEOMETRY_TYPES = dict( # Geometry related functions and classes follow. -cdef void * _createOgrGeomFromWKB(object wkb) except NULL: - """Make an OGR geometry from a WKB string""" - - geom_type = bytearray(wkb)[1] - cdef unsigned char *buffer = wkb - cdef void *cogr_geometry = _ogr.OGR_G_CreateGeometry(geom_type) - if cogr_geometry != NULL: - _ogr.OGR_G_ImportFromWkb(cogr_geometry, buffer, len(wkb)) - return cogr_geometry - cdef _deleteOgrGeom(void *cogr_geometry): """Delete an OGR geometry""" @@ -445,15 +470,6 @@ cdef class GeomBuilder: self.geom = geom return getattr(self, '_build' + self.geomtypename)() - cpdef build_wkb(self, object wkb): - """Builds a GeoJSON object from a Well-Known Binary format (WKB).""" - # The only other method anyone needs to call - cdef object data = wkb - cdef void *cogr_geometry = _createOgrGeomFromWKB(data) - result = self.build(cogr_geometry) - _deleteOgrGeom(cogr_geometry) - return result - cdef class OGRGeomBuilder: """ diff --git a/rasterio/dtypes.py b/rasterio/dtypes.py index 90097d81..8b2de802 100644 --- a/rasterio/dtypes.py +++ b/rasterio/dtypes.py @@ -74,6 +74,7 @@ def _gdal_typename(dt): except KeyError: return typename_fwd[dtype_rev[dt().dtype.name]] + def check_dtype(dt): if dt not in dtype_rev: try: @@ -83,32 +84,105 @@ def check_dtype(dt): return True -def get_minimum_int_dtype(values): +def get_minimum_dtype(values): """ - Uses range checking to determine the minimum integer data type required + Uses range checking to determine the minimum integer or floating point + data type required to represent values. - :param values: numpy array - :return: named data type that can be later used to create a numpy dtype + Parameters + ---------- + values: list-like + + + Returns + ------- + rasterio dtype string """ + import numpy + + if not is_ndarray(values): + values = numpy.array(values) + min_value = values.min() max_value = values.max() - - if min_value >= 0: - if max_value <= 255: - return uint8 - elif max_value <= 65535: - return uint16 - elif max_value <= 4294967295: - return uint32 - elif min_value >= -32768 and max_value <= 32767: - return int16 - elif min_value >= -2147483648 and max_value <= 2147483647: - return int32 + + if values.dtype.kind == 'i': + if min_value >= 0: + if max_value <= 255: + return uint8 + elif max_value <= 65535: + return uint16 + elif max_value <= 4294967295: + return uint32 + elif min_value >= -32768 and max_value <= 32767: + return int16 + elif min_value >= -2147483648 and max_value <= 2147483647: + return int32 + + else: + if min_value >= -3.4028235e+38 and max_value <= 3.4028235e+38: + return float32 + return float64 def is_ndarray(array): import numpy return isinstance(array, numpy.ndarray) or hasattr(array, '__array__') + + +def can_cast_dtype(values, dtype): + """ + Tests if values can be cast to dtype without loss of information. + + Parameters + ---------- + values: list-like + dtype: numpy dtype or string + + Returns + ------- + boolean + True if values can be cast to data type. + """ + + import numpy + + if not is_ndarray(values): + values = numpy.array(values) + + if values.dtype.name == numpy.dtype(dtype).name: + return True + + elif values.dtype.kind == 'f': + return numpy.allclose(values, values.astype(dtype)) + + else: + return numpy.array_equal(values, values.astype(dtype)) + + +def validate_dtype(values, valid_dtypes): + """ + Tests if dtype of values is one of valid_dtypes. + + Parameters + ---------- + values: list-like + valid_dtypes: list-like + list of valid dtype strings, e.g., ('int16', 'int32') + + Returns + ------- + boolean: + True if dtype of values is one of valid_dtypes + """ + + import numpy + + if not is_ndarray(values): + values = numpy.array(values) + + return (values.dtype.name in valid_dtypes or + get_minimum_dtype(values) in valid_dtypes) \ No newline at end of file diff --git a/rasterio/features.py b/rasterio/features.py index d52da1a2..d2aa7974 100644 --- a/rasterio/features.py +++ b/rasterio/features.py @@ -10,7 +10,7 @@ import numpy as np import rasterio from rasterio._features import _shapes, _sieve, _rasterize, _bounds from rasterio.transform import IDENTITY, guard_transform -from rasterio.dtypes import get_minimum_int_dtype +from rasterio.dtypes import validate_dtype, can_cast_dtype, get_minimum_dtype log = logging.getLogger('rasterio') @@ -104,18 +104,6 @@ def shapes(image, mask=None, connectivity=4, transform=IDENTITY): """ - valid_dtypes = ('int16', 'int32', 'uint8', 'uint16', 'float32') - - if np.dtype(image.dtype).name not in valid_dtypes: - raise ValueError('image dtype must be one of: %s' - % (', '.join(valid_dtypes))) - - if mask is not None and np.dtype(mask.dtype).name not in ('bool', 'uint8'): - raise ValueError("Mask must be dtype rasterio.bool_ or rasterio.uint8") - - if connectivity not in (4, 8): - raise ValueError("Connectivity Option must be 4 or 8") - transform = guard_transform(transform) with rasterio.drivers(): @@ -165,49 +153,18 @@ def sieve(image, size, out=None, output=None, mask=None, connectivity=4): """ - valid_dtypes = ('int16', 'int32', 'uint8', 'uint16') - - if np.dtype(image.dtype).name not in valid_dtypes: - valid_types_str = ', '.join(('rasterio.{0}'.format(t) for t - in valid_dtypes)) - raise ValueError('image dtype must be one of: %s' % valid_types_str) - - if size <= 0: - raise ValueError('size must be greater than 0') - elif type(size) == float: - raise ValueError('size must be an integer number of pixels') - elif size > (image.shape[0] * image.shape[1]): - raise ValueError('size must be smaller than size of image') - - if connectivity not in (4, 8): - raise ValueError('connectivity must be 4 or 8') - - if mask is not None: - if np.dtype(mask.dtype) not in ('bool', 'uint8'): - raise ValueError('Mask must be dtype rasterio.bool_ or ' - 'rasterio.uint8') - elif mask.shape != image.shape: - raise ValueError('mask shape must be same as image shape') - # Start moving users over to 'out'. if output is not None: warnings.warn( "The 'output' keyword arg has been superceded by 'out' " "and will be removed before Rasterio 1.0.", FutureWarning, - stacklevel=2) + stacklevel=2) # pragma: no cover out = out if out is not None else output + if out is None: - if isinstance(image, tuple): - out = np.zeros(image.shape, image.dtype) - else: - out = np.zeros_like(image) - else: - if np.dtype(image.dtype).name != np.dtype(out.dtype).name: - raise ValueError('out raster must match dtype of image') - elif out.shape != image.shape: - raise ValueError('out raster shape must be same as image shape') + out = np.zeros(image.shape, image.dtype) with rasterio.drivers(): _sieve(image, size, out, mask, connectivity) @@ -268,97 +225,92 @@ def rasterize( """ - valid_dtypes = ('int16', 'int32', 'uint8', 'uint16', 'uint32', 'float32', - 'float64') + valid_dtypes = ( + 'int16', 'int32', 'uint8', 'uint16', 'uint32', 'float32', 'float64' + ) - def get_valid_dtype(values): - values_dtype = values.dtype - if values_dtype.kind == 'i': - values_dtype = np.dtype(get_minimum_int_dtype(values)) - if values_dtype.name in valid_dtypes: - return values_dtype - return None + def format_invalid_dtype(param): + return '{0} dtype must be one of: {1}'.format( + param, ', '.join(valid_dtypes) + ) + + def format_cast_error(param, dtype): + return '{0} cannot be cast to specified dtype: {1}'.format(param, dtype) - def can_cast_dtype(values, dtype): - if values.dtype.name == np.dtype(dtype).name: - return True - elif values.dtype.kind == 'f': - return np.allclose(values, values.astype(dtype)) - else: - return np.array_equal(values, values.astype(dtype)) if fill != 0: fill_array = np.array([fill]) - if get_valid_dtype(fill_array) is None: - raise ValueError('fill must be one of these types: %s' - % (', '.join(valid_dtypes))) - elif dtype is not None and not can_cast_dtype(fill_array, dtype): - raise ValueError('fill value cannot be cast to specified dtype') + if not validate_dtype(fill_array, valid_dtypes): + raise ValueError(format_invalid_dtype('fill')) + + if dtype is not None and not can_cast_dtype(fill_array, dtype): + raise ValueError(format_cast_error('fill', dtype)) if default_value != 1: default_value_array = np.array([default_value]) - if get_valid_dtype(default_value_array) is None: - raise ValueError('default_value must be one of these types: %s' - % (', '.join(valid_dtypes))) - elif dtype is not None and not can_cast_dtype(default_value_array, - dtype): - raise ValueError('default_value cannot be cast to specified dtype') + if not validate_dtype(default_value_array, valid_dtypes): + raise ValueError(format_invalid_dtype('default_value')) + + if dtype is not None and not can_cast_dtype(default_value_array, dtype): + raise ValueError(format_cast_error('default_vaue', dtype)) + + if dtype is not None and np.dtype(dtype).name not in valid_dtypes: + raise ValueError(format_invalid_dtype('dtype')) + valid_shapes = [] shape_values = [] for index, item in enumerate(shapes): - try: - if isinstance(item, (tuple, list)): - geom, value = item - else: - geom = item - value = default_value - geom = getattr(geom, '__geo_interface__', None) or geom - if (not isinstance(geom, dict) or - 'type' not in geom or 'coordinates' not in geom): - raise ValueError( - 'Object %r at index %d is not a geometry object' % - (geom, index)) + if isinstance(item, (tuple, list)): + geom, value = item + else: + geom = item + value = default_value + geom = getattr(geom, '__geo_interface__', None) or geom + + #not isinstance(geom, dict) or + if 'type' in geom or 'coordinates' in geom: valid_shapes.append((geom, value)) shape_values.append(value) - except Exception: - log.exception('Exception caught, skipping shape %d', index) + + else: + raise ValueError( + 'Invalid geometry object at index {0}'.format(index) + ) if not valid_shapes: - raise ValueError('No valid shapes found for rasterize. Shapes must be ' - 'valid geometry objects') + raise ValueError('No valid geometry objects found for rasterize') shape_values = np.array(shape_values) - values_dtype = get_valid_dtype(shape_values) - if values_dtype is None: - raise ValueError('shape values must be one of these dtypes: %s' % - (', '.join(valid_dtypes))) + + if not validate_dtype(shape_values, valid_dtypes): + raise ValueError(format_invalid_dtype('shape values')) if dtype is None: - dtype = values_dtype - elif np.dtype(dtype).name not in valid_dtypes: - raise ValueError('dtype must be one of: %s' % (', '.join(valid_dtypes))) + dtype = get_minimum_dtype(np.append(shape_values, fill)) + elif not can_cast_dtype(shape_values, dtype): - raise ValueError('shape values could not be cast to specified dtype') + raise ValueError(format_cast_error('shape values', dtype)) if output is not None: warnings.warn( "The 'output' keyword arg has been superceded by 'out' " "and will be removed before Rasterio 1.0.", FutureWarning, - stacklevel=2) + stacklevel=2) # pragma: no cover + out = out if out is not None else output if out is not None: if np.dtype(out.dtype).name not in valid_dtypes: - raise ValueError('Output image dtype must be one of: %s' - % (', '.join(valid_dtypes))) + raise ValueError(format_invalid_dtype('out')) + if not can_cast_dtype(shape_values, out.dtype): - raise ValueError('shape values cannot be cast to dtype of output ' - 'image') + raise ValueError(format_cast_error('shape values', out.dtype.name)) elif out_shape is not None: out = np.empty(out_shape, dtype=dtype) out.fill(fill) + else: raise ValueError('Either an output shape or image must be provided') diff --git a/rasterio/rio/features.py b/rasterio/rio/features.py index 0345bbbd..dee16b23 100644 --- a/rasterio/rio/features.py +++ b/rasterio/rio/features.py @@ -122,14 +122,26 @@ def mask( bounds = geojson.get('bbox', calculate_bounds(geojson)) with rasterio.open(input) as src: - has_disjoint_bounds = disjoint_bounds(bounds, src.bounds) + # If y pixel value is positive, then invert y dimension in bounds + invert_y = src.affine.e > 0 + + src_bounds = src.bounds + if invert_y: + src_bounds = [src.bounds[0], src.bounds[3], + src.bounds[2], src.bounds[1]] + + has_disjoint_bounds = disjoint_bounds(bounds, src_bounds) if crop: if has_disjoint_bounds: + raise click.BadParameter('not allowed for GeoJSON outside ' 'the extent of the input raster', param=crop, param_hint='--crop') + if invert_y: + bounds = (bounds[0], bounds[3], bounds[2], bounds[1]) + window = src.window(*bounds) transform = src.window_transform(window) (r1, r2), (c1, c2) = window @@ -256,6 +268,9 @@ def shapes( def __call__(self): with rasterio.open(input) as src: + if bidx is not None and bidx > src.count: + raise ValueError('bidx is out of range for raster') + img = None msk = None @@ -533,6 +548,12 @@ def rasterize( kwargs = template_ds.meta.copy() kwargs['count'] = 1 + + # DEPRECATED + # upgrade transform to affine object or we may get an invalid + # transform set on output + kwargs['transform'] = template_ds.affine + template_ds.close() else: diff --git a/tests/conftest.py b/tests/conftest.py index 90aed75a..16f39e73 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -7,6 +7,10 @@ import sys from click.testing import CliRunner import py import pytest +import numpy + + +DEFAULT_SHAPE = (10, 10) if sys.version_info > (3,): @@ -38,3 +42,163 @@ def data(): for filename in test_files: shutil.copy(filename, str(tmpdir)) return tmpdir + + +@pytest.fixture +def basic_geometry(): + """ + Returns + ------- + + dict: GeoJSON-style geometry object. + Coordinates are in grid coordinates (Affine.identity()). + """ + + return { + 'type': 'Polygon', + 'coordinates': [[(2, 2), (2, 4.25), (4.25, 4.25), (4.25, 2), (2, 2)]] + } + + +@pytest.fixture +def basic_feature(basic_geometry): + """ + Returns + ------- + + dict: GeoJSON object. + Coordinates are in grid coordinates (Affine.identity()). + """ + + return { + 'geometry': basic_geometry, + 'properties': { + 'val': 15 + }, + 'type': 'Feature' + } + + +@pytest.fixture +def basic_featurecollection(basic_feature): + """ + Returns + ------- + + dict: GeoJSON FeatureCollection object. + Coordinates are in grid coordinates (Affine.identity()). + """ + + return { + 'features': [basic_feature], + 'type': 'FeatureCollection' + } + + +@pytest.fixture +def basic_image(): + """ + A basic 10x10 array for testing sieve and shapes functions. + Contains a square feature 3x3 (size 9). + Equivalent to results of rasterizing basic_geometry with all_touched=True. + + Returns + ------- + + numpy ndarray + """ + + image = numpy.zeros(DEFAULT_SHAPE, dtype=numpy.uint8) + image[2:5, 2:5] = 1 + + return image + + +@pytest.fixture +def basic_image_2x2(): + """ + A basic 10x10 array for testing sieve and shapes functions. + Contains a square feature 2x2 (size 4). + Equivalent to results of rasterizing basic_geometry with all_touched=False. + + Returns + ------- + + numpy ndarray + """ + + image = numpy.zeros(DEFAULT_SHAPE, dtype=numpy.uint8) + image[2:4, 2:4] = 1 + + return image + + +@pytest.fixture +def pixelated_image(basic_image): + """ + A basic 10x10 array for testing sieve functions. Contains a square feature + 3x3 (size 9), with 2 isolated pixels. + + Returns + ------- + + numpy ndarray + """ + + image = basic_image.copy() + image [0, 0] = 1 + image [8, 8] = 1 + + return image + + +@pytest.fixture +def diagonal_image(): + """ + A 10x10 array for testing sieve functions, with only one diagonal filled. + + Returns + ------- + + numpy ndarray + """ + + image = numpy.zeros(DEFAULT_SHAPE, dtype=numpy.uint8) + numpy.fill_diagonal(image, 1) + return image + + +@pytest.fixture() +def pixelated_image_file(tmpdir, pixelated_image): + """ + A basic raster file with a 10x10 array for testing sieve functions. + Contains data from pixelated_image. + + Returns + ------- + + string + Filename of test raster file + """ + + from affine import Affine + import rasterio + + image = pixelated_image + + outfilename = str(tmpdir.join('pixelated_image.tif')) + kwargs = { + "crs": {'init': 'epsg:4326'}, + "transform": Affine.identity(), + "count": 1, + "dtype": rasterio.uint8, + "driver": "GTiff", + "width": image.shape[1], + "height": image.shape[0], + "nodata": 255 + } + with rasterio.drivers(): + with rasterio.open(outfilename, 'w', **kwargs) as out: + out.write_band(1, image) + + return outfilename diff --git a/tests/test_dtypes.py b/tests/test_dtypes.py index 1826ec52..9c0310b0 100644 --- a/tests/test_dtypes.py +++ b/tests/test_dtypes.py @@ -1,23 +1,58 @@ import numpy as np -from rasterio import dtypes, ubyte +from rasterio import ( + ubyte, uint8, uint16, uint32, int16, int32, float32, float64) +from rasterio.dtypes import ( + _gdal_typename, is_ndarray, check_dtype, get_minimum_dtype, can_cast_dtype, + validate_dtype +) def test_is_ndarray(): - assert dtypes.is_ndarray(np.zeros((1,))) - assert dtypes.is_ndarray([0]) == False - assert dtypes.is_ndarray((0,)) == False + assert is_ndarray(np.zeros((1,))) + assert is_ndarray([0]) == False + assert is_ndarray((0,)) == False def test_np_dt_uint8(): - assert dtypes.check_dtype(np.uint8) + assert check_dtype(np.uint8) def test_dt_ubyte(): - assert dtypes.check_dtype(ubyte) + assert check_dtype(ubyte) + + +def test_check_dtype_invalid(): + assert check_dtype('foo') == False def test_gdal_name(): - assert dtypes._gdal_typename(ubyte) == 'Byte' - assert dtypes._gdal_typename(np.uint8) == 'Byte' - assert dtypes._gdal_typename(np.uint16) == 'UInt16' + assert _gdal_typename(ubyte) == 'Byte' + assert _gdal_typename(np.uint8) == 'Byte' + assert _gdal_typename(np.uint16) == 'UInt16' + + +def test_get_minimum_dtype(): + assert get_minimum_dtype([0, 1]) == uint8 + assert get_minimum_dtype([0, 1000]) == uint16 + assert get_minimum_dtype([0, 100000]) == uint32 + assert get_minimum_dtype([-1, 0, 1]) == int16 + assert get_minimum_dtype([-1, 0, 100000]) == int32 + assert get_minimum_dtype([-1.5, 0, 1.5]) == float32 + assert get_minimum_dtype([-1.5e+100, 0, 1.5e+100]) == float64 + + +def test_can_cast_dtype(): + assert can_cast_dtype((1, 2, 3), np.uint8) == True + assert can_cast_dtype(np.array([1, 2, 3]), np.uint8) == True + assert can_cast_dtype(np.array([1, 2, 3], dtype=np.uint8), np.uint8) == True + assert can_cast_dtype(np.array([1, 2, 3]), np.float32) == True + assert can_cast_dtype(np.array([1.4, 2.1, 3.65]), np.float32) == True + assert can_cast_dtype(np.array([1.4, 2.1, 3.65]), np.uint8) == False + + +def test_validate_dtype(): + assert validate_dtype([1, 2, 3], ('uint8', 'uint16')) == True + assert validate_dtype(np.array([1, 2, 3]), ('uint8', 'uint16')) == True + assert validate_dtype(np.array([1.4, 2.1, 3.65]), ('float32',)) == True + assert validate_dtype(np.array([1.4, 2.1, 3.65]),('uint8',)) == False diff --git a/tests/test_features.py b/tests/test_features.py new file mode 100644 index 00000000..bdbdab22 --- /dev/null +++ b/tests/test_features.py @@ -0,0 +1,728 @@ +import logging +import sys +import numpy +import pytest + +from affine import Affine +import rasterio +from rasterio.features import bounds, geometry_mask, rasterize, sieve, shapes + + +DEFAULT_SHAPE = (10, 10) + +logging.basicConfig(stream=sys.stderr, level=logging.DEBUG) + + +def test_bounds_point(): + g = {'type': 'Point', 'coordinates': [10, 10]} + assert bounds(g) == (10, 10, 10, 10) + + +def test_bounds_line(): + g = {'type': 'LineString', 'coordinates': [[0, 0], [10, 10]]} + assert bounds(g) == (0, 0, 10, 10) + + +def test_bounds_polygon(): + g = {'type': 'Polygon', 'coordinates': [[[0, 0], [10, 10], [10, 0]]]} + assert bounds(g) == (0, 0, 10, 10) + + +def test_bounds_z(): + g = {'type': 'Point', 'coordinates': [10, 10, 10]} + assert bounds(g) == (10, 10, 10, 10) + + +def test_bounds_invalid_obj(): + with pytest.raises(KeyError): + bounds({'type': 'bogus', 'not_coordinates': []}) + + +def test_feature_collection(basic_featurecollection): + fc = basic_featurecollection + assert bounds(fc) == bounds(fc['features'][0]) == (2, 2, 4.25, 4.25) + + +def test_bounds_existing_bbox(basic_featurecollection): + """ + Test with existing bbox in geojson, similar to that produced by + rasterio. Values specifically modified here for testing, bboxes are not + valid as written. + """ + + fc = basic_featurecollection + fc['bbox'] = [0, 10, 10, 20] + fc['features'][0]['bbox'] = [0, 100, 10, 200] + + assert bounds(fc['features'][0]) == (0, 100, 10, 200) + assert bounds(fc) == (0, 10, 10, 20) + + +def test_geometry_mask(basic_geometry, basic_image_2x2): + with rasterio.drivers(): + assert numpy.array_equal( + basic_image_2x2 == 0, + geometry_mask( + [basic_geometry], + out_shape=DEFAULT_SHAPE, + transform=Affine.identity() + ) + ) + + +def test_geometry_mask_invert(basic_geometry, basic_image_2x2): + with rasterio.drivers(): + assert numpy.array_equal( + basic_image_2x2, + geometry_mask( + [basic_geometry], + out_shape=DEFAULT_SHAPE, + transform=Affine.identity(), + invert=True + ) + ) + + +def test_rasterize(basic_geometry, basic_image_2x2): + """ Rasterize operation should succeed for both an out_shape and out """ + + with rasterio.drivers(): + assert numpy.array_equal( + basic_image_2x2, + rasterize([basic_geometry], out_shape=DEFAULT_SHAPE) + ) + + out = numpy.zeros(DEFAULT_SHAPE) + rasterize([basic_geometry], out=out) + assert numpy.array_equal(basic_image_2x2, out) + + +def test_rasterize_invalid_out_dtype(basic_geometry): + """ A non-supported data type for out should raise an exception """ + + out = numpy.zeros(DEFAULT_SHAPE, dtype=numpy.int64) + with rasterio.drivers(): + with pytest.raises(ValueError): + rasterize([basic_geometry], out=out) + + +def test_rasterize_shapes_out_dtype_mismatch(basic_geometry): + """ Shape values must be able to fit in data type for out """ + + out = numpy.zeros(DEFAULT_SHAPE, dtype=numpy.uint8) + with rasterio.drivers(): + with pytest.raises(ValueError): + rasterize([(basic_geometry, 10000000)], out=out) + + +def test_rasterize_missing_out(basic_geometry): + """ If both out and out_shape are missing, should raise exception """ + + with rasterio.drivers(): + with pytest.raises(ValueError): + rasterize([basic_geometry], out=None, out_shape=None) + + +def test_rasterize_missing_shapes(): + """ Shapes are required for this operation """ + + with rasterio.drivers(): + with pytest.raises(ValueError) as ex: + rasterize([], out_shape=DEFAULT_SHAPE) + + assert 'No valid geometry objects' in str(ex.value) + + +def test_rasterize_invalid_shapes(): + """ Invalid shapes should raise an exception rather than be skipped """ + + with rasterio.drivers(): + with pytest.raises(ValueError) as ex: + rasterize([{'foo': 'bar'}], out_shape=DEFAULT_SHAPE) + + assert 'Invalid geometry object' in str(ex.value) + + +def test_rasterize_default_value(basic_geometry, basic_image_2x2): + """ All shapes should rasterize to the default value """ + + default_value = 2 + truth = basic_image_2x2 * default_value + + with rasterio.drivers(): + assert numpy.array_equal( + truth, + rasterize( + [basic_geometry], out_shape=DEFAULT_SHAPE, + default_value=default_value + ) + ) + + +def test_rasterize_invalid_default_value(basic_geometry): + """ A default value that requires an int64 should raise an exception """ + + with rasterio.drivers(): + with pytest.raises(ValueError): + rasterize( + [basic_geometry], out_shape=DEFAULT_SHAPE, + default_value=1000000000000 + ) + + +def test_rasterize_fill_value(basic_geometry, basic_image_2x2): + """ All pixels not covered by shapes should be given fill value """ + + default_value = 2 + with rasterio.drivers(): + assert numpy.array_equal( + basic_image_2x2 + 1, + rasterize( + [basic_geometry], out_shape=DEFAULT_SHAPE, fill=1, + default_value=default_value + ) + ) + + +def test_rasterize_invalid_fill_value(basic_geometry): + """ A fill value that requires an int64 should raise an exception """ + + with rasterio.drivers(): + with pytest.raises(ValueError): + rasterize( + [basic_geometry], out_shape=DEFAULT_SHAPE, fill=1000000000000, + default_value=2 + ) + + +def test_rasterize_fill_value_dtype_mismatch(basic_geometry): + """ A fill value that doesn't match dtype should fail """ + + with rasterio.drivers(): + with pytest.raises(ValueError): + rasterize( + [basic_geometry], out_shape=DEFAULT_SHAPE, fill=1000000, + default_value=2, dtype=numpy.uint8 + ) + + +def test_rasterize_all_touched(basic_geometry, basic_image): + with rasterio.drivers(): + assert numpy.array_equal( + basic_image, + rasterize( + [basic_geometry], out_shape=DEFAULT_SHAPE, all_touched=True + ) + ) + + +def test_rasterize_value(basic_geometry, basic_image_2x2): + """ + All shapes should rasterize to the value passed in a tuple alongside + each shape + """ + + value = 5 + with rasterio.drivers(): + assert numpy.array_equal( + basic_image_2x2 * value, + rasterize( + [(basic_geometry, value)], out_shape=DEFAULT_SHAPE + ) + ) + + +def test_rasterize_invalid_value(basic_geometry): + """ A shape value that requires an int64 should raise an exception """ + + with rasterio.drivers(): + with pytest.raises(ValueError) as ex: + rasterize( + [(basic_geometry, 1000000000000)], out_shape=DEFAULT_SHAPE + ) + + assert 'dtype must be one of' in str(ex.value) + + +def test_rasterize_supported_dtype(basic_geometry): + """ Supported data types should return valid results """ + + with rasterio.drivers(): + supported_types = ( + ('int16', -32768), + ('int32', -2147483648), + ('uint8', 255), + ('uint16', 65535), + ('uint32', 4294967295), + ('float32', 1.434532), + ('float64', -98332.133422114) + ) + + for dtype, default_value in supported_types: + truth = numpy.zeros(DEFAULT_SHAPE, dtype=dtype) + truth[2:4, 2:4] = default_value + + result = rasterize( + [basic_geometry], + out_shape=DEFAULT_SHAPE, + default_value=default_value, + dtype=dtype + ) + assert numpy.array_equal(result, truth) + assert numpy.dtype(result.dtype) == numpy.dtype(truth.dtype) + + result = rasterize( + [(basic_geometry, default_value)], + out_shape=DEFAULT_SHAPE + ) + if numpy.dtype(dtype).kind == 'f': + assert numpy.allclose(result, truth) + else: + assert numpy.array_equal(result, truth) + # Since dtype is auto-detected, it may not match due to upcasting + + +def test_rasterize_unsupported_dtype(basic_geometry): + """ Unsupported types should all raise exceptions """ + + with rasterio.drivers(): + unsupported_types = ( + ('int8', -127), + ('int64', 20439845334323), + ('float16', -9343.232) + ) + + for dtype, default_value in unsupported_types: + with pytest.raises(ValueError): + rasterize( + [basic_geometry], + out_shape=DEFAULT_SHAPE, + default_value=default_value, + dtype=dtype + ) + + with pytest.raises(ValueError): + rasterize( + [(basic_geometry, default_value)], + out_shape=DEFAULT_SHAPE, + dtype=dtype + ) + + +def test_rasterize_mismatched_dtype(basic_geometry): + """ Mismatched values and dtypes should raise exceptions """ + + with rasterio.drivers(): + mismatched_types = (('uint8', 3.2423), ('uint8', -2147483648)) + for dtype, default_value in mismatched_types: + with pytest.raises(ValueError): + rasterize( + [basic_geometry], + out_shape=DEFAULT_SHAPE, + default_value=default_value, + dtype=dtype + ) + + with pytest.raises(ValueError): + rasterize( + [(basic_geometry, default_value)], + out_shape=DEFAULT_SHAPE, + dtype=dtype + ) + + +def test_rasterize_geometries_symmetric(): + """ Make sure that rasterize is symmetric with shapes """ + + transform = (1.0, 0.0, 0.0, 0.0, -1.0, 0.0) + truth = numpy.zeros(DEFAULT_SHAPE, dtype=rasterio.ubyte) + truth[2:5, 2:5] = 1 + with rasterio.drivers(): + s = shapes(truth, transform=transform) + result = rasterize(s, out_shape=DEFAULT_SHAPE, transform=transform) + assert numpy.array_equal(result, truth) + + +def test_rasterize_internal_driver_manager(basic_geometry): + """ Rasterize should work without explicitly calling driver manager """ + + assert rasterize([basic_geometry], out_shape=DEFAULT_SHAPE).sum() == 4 + + +def test_shapes(basic_image): + """ Test creation of shapes from pixel values """ + + with rasterio.drivers(): + results = list(shapes(basic_image)) + + assert len(results) == 2 + + shape, value = results[0] + assert shape == { + 'coordinates': [ + [(2, 2), (2, 5), (5, 5), (5, 2), (2, 2)] + ], + 'type': 'Polygon' + } + assert value == 1 + + shape, value = results[1] + assert shape == { + 'coordinates': [ + [(0, 0), (0, 10), (10, 10), (10, 0), (0, 0)], + [(2, 2), (5, 2), (5, 5), (2, 5), (2, 2)] + ], + 'type': 'Polygon' + } + assert value == 0 + + +def test_shapes_band(pixelated_image, pixelated_image_file): + """ Shapes from a band should match shapes from an array """ + + with rasterio.drivers(): + truth = list(shapes(pixelated_image)) + + with rasterio.open(pixelated_image_file) as src: + band = rasterio.band(src, 1) + assert truth == list(shapes(band)) + + # Mask band should function, but will mask out some results + assert truth[0] == list(shapes(band, mask=band))[0] + + +def test_shapes_connectivity_rook(diagonal_image): + """ + Diagonals are not connected, so there will be 1 feature per pixel plus + background. + """ + + with rasterio.drivers(): + assert len(list(shapes(diagonal_image, connectivity=4))) == 12 + + +def test_shapes_connectivity_queen(diagonal_image): + """ + Diagonals are connected, so there will be 1 feature for all pixels plus + background. + """ + + with rasterio.drivers(): + assert len(list(shapes(diagonal_image, connectivity=8))) == 2 + + +def test_shapes_connectivity_invalid(diagonal_image): + """ Invalid connectivity should raise exception """ + + with rasterio.drivers(): + with pytest.raises(ValueError): + assert next(shapes(diagonal_image, connectivity=12)) + + +def test_shapes_mask(basic_image): + """ Only pixels not masked out should be converted to features """ + + mask = numpy.ones(basic_image.shape, dtype=rasterio.bool_) + mask[4:5, 4:5] = False + + with rasterio.drivers(): + results = list(shapes(basic_image, mask=mask)) + + assert len(results) == 2 + + shape, value = results[0] + assert shape == { + 'coordinates': [ + [(2, 2), (2, 5), (4, 5), (4, 4), (5, 4), (5, 2), (2, 2)] + ], + 'type': 'Polygon' + } + assert value == 1 + + +def test_shapes_blank_mask(basic_image): + """ Mask is blank so results should mask shapes without mask """ + + with rasterio.drivers(): + assert numpy.array_equal( + list(shapes( + basic_image, + mask=numpy.ones(basic_image.shape, dtype=rasterio.bool_)) + ), + list(shapes(basic_image)) + ) + + +def test_shapes_invalid_mask_shape(basic_image): + """ A mask that is the wrong shape should fail """ + + with rasterio.drivers(): + with pytest.raises(ValueError): + next(shapes( + basic_image, + mask=numpy.ones( + (basic_image.shape[0] + 10, basic_image.shape[1] + 10), + dtype=rasterio.bool_ + ) + )) + + +def test_shapes_invalid_mask_dtype(basic_image): + """ A mask that is the wrong dtype should fail """ + + with rasterio.drivers(): + for dtype in ('int8', 'int16', 'int32'): + with pytest.raises(ValueError): + next(shapes( + basic_image, + mask=numpy.ones(basic_image.shape, dtype=dtype) + )) + + +def test_shapes_supported_dtypes(basic_image): + """ Supported data types should return valid results """ + + supported_types = ( + ('int16', -32768), + ('int32', -2147483648), + ('uint8', 255), + ('uint16', 65535), + ('float32', 1.434532) + ) + + with rasterio.drivers(): + for dtype, test_value in supported_types: + shape, value = next(shapes(basic_image.astype(dtype) * test_value)) + assert numpy.allclose(value, test_value) + + +def test_shapes_unsupported_dtypes(basic_image): + """ Unsupported data types should raise exceptions """ + + unsupported_types = ( + ('int8', -127), + ('uint32', 4294967295), + ('int64', 20439845334323), + ('float16', -9343.232), + ('float64', -98332.133422114) + ) + + with rasterio.drivers(): + for dtype, test_value in unsupported_types: + with pytest.raises(ValueError): + next(shapes(basic_image.astype(dtype) * test_value)) + + +def test_shapes_internal_driver_manager(basic_image): + """ Shapes should work without explicitly calling driver manager """ + + assert next(shapes(basic_image))[0]['type'] == 'Polygon' + + +def test_sieve_small(basic_image, pixelated_image): + """ + Setting the size smaller than or equal to the size of the feature in the + image should not change the image. + """ + + with rasterio.drivers(): + assert numpy.array_equal( + basic_image, + sieve(pixelated_image, basic_image.sum()) + ) + + +def test_sieve_large(basic_image): + """ + Setting the size larger than size of feature should leave us an empty image. + """ + + with rasterio.drivers(): + assert not numpy.any(sieve(basic_image, basic_image.sum() + 1)) + + +def test_sieve_invalid_size(basic_image): + with rasterio.drivers(): + for invalid_size in (0, 45.1234, basic_image.size + 1): + with pytest.raises(ValueError): + sieve(basic_image, invalid_size) + + +def test_sieve_connectivity_rook(diagonal_image): + """ Diagonals are not connected, so feature is removed """ + + assert not numpy.any( + sieve(diagonal_image, diagonal_image.sum(), connectivity=4) + ) + + +def test_sieve_connectivity_queen(diagonal_image): + """ Diagonals are connected, so feature is retained """ + + assert numpy.array_equal( + diagonal_image, + sieve(diagonal_image, diagonal_image.sum(), connectivity=8) + ) + + +def test_sieve_connectivity_invalid(basic_image): + with pytest.raises(ValueError): + sieve(basic_image, 54, connectivity=12) + + +def test_sieve_out(basic_image): + """ Output array passed in should match the returned array """ + + with rasterio.drivers(): + output = numpy.zeros_like(basic_image) + output[1:3, 1:3] = 5 + sieved_image = sieve(basic_image, basic_image.sum(), out=output) + assert numpy.array_equal(basic_image, sieved_image) + assert numpy.array_equal(output, sieved_image) + + +def test_sieve_invalid_out(basic_image): + """ Output with different dtype or shape should fail """ + + with rasterio.drivers(): + with pytest.raises(ValueError): + sieve( + basic_image, basic_image.sum(), + out=numpy.zeros(basic_image.shape, dtype=rasterio.int32) + ) + + with pytest.raises(ValueError): + sieve( + basic_image, basic_image.sum(), + out=numpy.zeros( + (basic_image.shape[0] + 10, basic_image.shape[1] + 10), + dtype=rasterio.ubyte + ) + ) + + +def test_sieve_mask(basic_image): + """ + Only areas within the overlap of mask and input will be kept, so long + as mask is a bool or uint8 dtype. + """ + + mask = numpy.ones(basic_image.shape, dtype=rasterio.bool_) + mask[4:5, 4:5] = False + truth = basic_image * numpy.invert(mask) + + with rasterio.drivers(): + sieved_image = sieve(basic_image, basic_image.sum(), mask=mask) + assert sieved_image.sum() > 0 + + assert numpy.array_equal( + truth, + sieved_image + ) + + assert numpy.array_equal( + truth.astype(rasterio.uint8), + sieved_image + ) + + +def test_sieve_blank_mask(basic_image): + """ A blank mask should have no effect """ + + mask = numpy.ones(basic_image.shape, dtype=rasterio.bool_) + with rasterio.drivers(): + assert numpy.array_equal( + basic_image, + sieve(basic_image, basic_image.sum(), mask=mask) + ) + + +def test_sieve_invalid_mask_shape(basic_image): + """ A mask that is the wrong shape should fail """ + + with rasterio.drivers(): + with pytest.raises(ValueError): + sieve( + basic_image, basic_image.sum(), + mask=numpy.ones( + (basic_image.shape[0] + 10, basic_image.shape[1] + 10), + dtype=rasterio.bool_ + ) + ) + + +def test_sieve_invalid_mask_dtype(basic_image): + """ A mask that is the wrong dtype should fail """ + + with rasterio.drivers(): + for dtype in ('int8', 'int16', 'int32'): + with pytest.raises(ValueError): + sieve( + basic_image, basic_image.sum(), + mask=numpy.ones(basic_image.shape, dtype=dtype) + ) + + +def test_sieve_supported_dtypes(basic_image): + """ Supported data types should return valid results """ + + supported_types = ( + ('int16', -32768), + ('int32', -2147483648), + ('uint8', 255), + ('uint16', 65535) + ) + + with rasterio.drivers(): + for dtype, test_value in supported_types: + truth = (basic_image).astype(dtype) * test_value + sieved_image = sieve(truth, basic_image.sum()) + assert numpy.array_equal(truth, sieved_image) + assert numpy.dtype(sieved_image.dtype) == numpy.dtype(dtype) + + +def test_sieve_unsupported_dtypes(basic_image): + """ Unsupported data types should raise exceptions """ + + unsupported_types = ( + ('int8', -127), + ('uint32', 4294967295), + ('int64', 20439845334323), + ('float16', -9343.232), + ('float32', 1.434532), + ('float64', -98332.133422114) + ) + + with rasterio.drivers(): + for dtype, test_value in unsupported_types: + with pytest.raises(ValueError): + sieve( + (basic_image).astype(dtype) * test_value, + basic_image.sum() + ) + + +def test_sieve_band(pixelated_image, pixelated_image_file): + """ Sieving a band from a raster file should match sieve of array """ + + with rasterio.drivers(): + truth = sieve(pixelated_image, 9) + + with rasterio.open(pixelated_image_file) as src: + band = rasterio.band(src, 1) + assert numpy.array_equal(truth, sieve(band, 9)) + + # Mask band should also work but will be a no-op + assert numpy.array_equal( + pixelated_image, + sieve(band, 9, mask=band) + ) + + +def test_sieve_internal_driver_manager(basic_image, pixelated_image): + """ Sieve should work without explicitly calling driver manager """ + + assert numpy.array_equal( + basic_image, + sieve(pixelated_image, basic_image.sum()) + ) diff --git a/tests/test_features_bounds.py b/tests/test_features_bounds.py deleted file mode 100644 index a7015cdf..00000000 --- a/tests/test_features_bounds.py +++ /dev/null @@ -1,65 +0,0 @@ -from rasterio.features import bounds - - -# Tests copied from Fiona 1.4.1 - -def test_bounds_point(): - g = {'type': 'Point', 'coordinates': [10, 10]} - assert bounds(g) == (10, 10, 10, 10) - - -def test_bounds_line(): - g = {'type': 'LineString', 'coordinates': [[0, 0], [10, 10]]} - assert bounds(g) == (0, 0, 10, 10) - - -def test_bounds_polygon(): - g = {'type': 'Polygon', 'coordinates': [[[0, 0], [10, 10], [10, 0]]]} - assert bounds(g) == (0, 0, 10, 10) - - -def test_bounds_z(): - g = {'type': 'Point', 'coordinates': [10, 10, 10]} - assert bounds(g) == (10, 10, 10, 10) - - -# TODO: add these to Fiona with update to bounds -def test_bounds_existing_bbox(): - """ Test with existing bbox in geojson, similar to that produced by - rasterio. Values specifically modified here for testing, bboxes are not - valid as written. - """ - - fc = { - 'bbox': [-107, 40, -105, 41], - 'features': [{ - 'bbox': [-107, 40, -104, 42], - 'geometry': { - 'coordinates': [ - [[-107, 40], [-106, 40], [-106, 41], [-107, 41], [-107, 40]] - ], - 'type': 'Polygon' - }, - 'type': 'Feature' - }], - 'type': 'FeatureCollection' - } - assert bounds(fc['features'][0]) == (-107, 40, -104, 42) - assert bounds(fc) == (-107, 40, -105, 41) - - -def test_feature_collection(): - fc = { - 'features': [{ - 'geometry': { - 'coordinates': [ - [[-107, 40], [-106, 40], [-106, 41], [-107, 41], [-107, 40]] - ], - 'type': 'Polygon' - }, - 'type': 'Feature' - }], - 'type': 'FeatureCollection' - } - assert bounds(fc['features'][0]) == (-107, 40, -106, 41) - assert bounds(fc) == (-107, 40, -106, 41) diff --git a/tests/test_features_rasterize.py b/tests/test_features_rasterize.py deleted file mode 100644 index 8848a8f5..00000000 --- a/tests/test_features_rasterize.py +++ /dev/null @@ -1,181 +0,0 @@ -import logging -import sys -import numpy -import pytest - -import rasterio -from rasterio.features import shapes, rasterize, geometry_mask - - -logging.basicConfig(stream=sys.stderr, level=logging.DEBUG) - - -def test_rasterize_geometries(): - """ - Make sure that geometries are correctly rasterized according to parameters - """ - - rows = cols = 10 - transform = (1.0, 0.0, 0.0, 0.0, 1.0, 0.0) - geometry = { - 'type': 'Polygon', - 'coordinates': [[(2, 2), (2, 4.25), (4.25, 4.25), (4.25, 2), (2, 2)]] - } - - with rasterio.drivers(): - # we expect a subset of the pixels using default mode - result = rasterize([geometry], out_shape=(rows, cols)) - truth = numpy.zeros((rows, cols)) - truth[2:4, 2:4] = 1 - assert numpy.array_equal(result, truth) - - out = numpy.zeros((rows, cols)) - result = rasterize([geometry], out=out, default_value=1) - assert numpy.array_equal(out, truth) - - # we expect all touched pixels - result = rasterize( - [geometry], out_shape=(rows, cols), all_touched=True - ) - truth = numpy.zeros((rows, cols)) - truth[2:5, 2:5] = 1 - assert numpy.array_equal(result, truth) - - # we expect the pixel value to match the one we pass in - value = 5 - result = rasterize([(geometry, value)], out_shape=(rows, cols)) - truth = numpy.zeros((rows, cols)) - truth[2:4, 2:4] = value - assert numpy.array_equal(result, truth) - - # Check the fill and default transform. - # we expect the pixel value to match the one we pass in - value = 5 - result = rasterize( - [(geometry, value)], - out_shape=(rows, cols), - fill=1 - ) - truth = numpy.ones((rows, cols)) - truth[2:4, 2:4] = value - assert numpy.array_equal(result, truth) - - -def test_rasterize_dtype(): - """Make sure that data types are handled correctly""" - - rows = cols = 10 - transform = (1.0, 0.0, 0.0, 0.0, 1.0, 0.0) - geometry = { - 'type': 'Polygon', - 'coordinates': [[(2, 2), (2, 4.25), (4.25, 4.25), (4.25, 2), (2, 2)]] - } - - with rasterio.drivers(): - # Supported types should all work properly - supported_types = ( - ('int16', -32768), - ('int32', -2147483648), - ('uint8', 255), - ('uint16', 65535), - ('uint32', 4294967295), - ('float32', 1.434532), - ('float64', -98332.133422114) - ) - - for dtype, default_value in supported_types: - truth = numpy.zeros((rows, cols), dtype=dtype) - truth[2:4, 2:4] = default_value - - result = rasterize( - [geometry], - out_shape=(rows, cols), - default_value=default_value, - dtype=dtype - ) - assert numpy.array_equal(result, truth) - assert numpy.dtype(result.dtype) == numpy.dtype(truth.dtype) - - result = rasterize( - [(geometry, default_value)], - out_shape=(rows, cols) - ) - if numpy.dtype(dtype).kind == 'f': - assert numpy.allclose(result, truth) - else: - assert numpy.array_equal(result, truth) - # Since dtype is auto-detected, it may not match due to upcasting - - # Unsupported types should all raise exceptions - unsupported_types = ( - ('int8', -127), - ('int64', 20439845334323), - ('float16', -9343.232) - ) - - for dtype, default_value in unsupported_types: - with pytest.raises(ValueError): - rasterize( - [geometry], - out_shape=(rows, cols), - default_value=default_value, - dtype=dtype - ) - - with pytest.raises(ValueError): - rasterize( - [(geometry, default_value)], - out_shape=(rows, cols), - dtype=dtype - ) - - # Mismatched values and dtypes should raise exceptions - mismatched_types = (('uint8', 3.2423), ('uint8', -2147483648)) - for dtype, default_value in mismatched_types: - with pytest.raises(ValueError): - rasterize( - [geometry], - out_shape=(rows, cols), - default_value=default_value, - dtype=dtype - ) - - with pytest.raises(ValueError): - rasterize( - [(geometry, default_value)], - out_shape=(rows, cols), - dtype=dtype - ) - - -def test_rasterize_geometries_symmetric(): - """Make sure that rasterize is symmetric with shapes""" - - rows = cols = 10 - transform = (1.0, 0.0, 0.0, 0.0, -1.0, 0.0) - truth = numpy.zeros((rows, cols), dtype=rasterio.ubyte) - truth[2:5, 2:5] = 1 - with rasterio.drivers(): - s = shapes(truth, transform=transform) - result = rasterize(s, out_shape=(rows, cols), transform=transform) - assert numpy.array_equal(result, truth) - - -def test_geometry_mask(): - rows = cols = 10 - transform = (1.0, 0.0, 0.0, 0.0, -1.0, 0.0) - truth = numpy.zeros((rows, cols), dtype=rasterio.bool_) - truth[2:5, 2:5] = True - with rasterio.drivers(): - s = shapes((truth * 10).astype(rasterio.ubyte), transform=transform) - # Strip out values returned from shapes, and only keep first shape - geoms = [next(s)[0]] - - # Regular mask should be the inverse of truth raster - mask = geometry_mask(geoms, out_shape=(rows, cols), transform=transform) - assert numpy.array_equal(mask, numpy.invert(truth)) - - # Inverted mask should be the same as the truth raster - mask = geometry_mask(geoms, out_shape=(rows, cols), transform=transform, - invert=True) - assert numpy.array_equal(mask, truth) \ No newline at end of file diff --git a/tests/test_features_shapes.py b/tests/test_features_shapes.py deleted file mode 100644 index 2419a150..00000000 --- a/tests/test_features_shapes.py +++ /dev/null @@ -1,144 +0,0 @@ -import logging -import sys -import numpy -import pytest - -import rasterio -import rasterio.features as ftrz - - -logging.basicConfig(stream=sys.stderr, level=logging.DEBUG) - - -def test_shapes(): - """Test creation of shapes from pixel values""" - - image = numpy.zeros((20, 20), dtype=rasterio.ubyte) - image[5:15, 5:15] = 127 - with rasterio.drivers(): - shapes = ftrz.shapes(image) - shape, val = next(shapes) - assert shape['type'] == 'Polygon' - assert len(shape['coordinates']) == 2 # exterior and hole - assert val == 0 - shape, val = next(shapes) - assert shape['type'] == 'Polygon' - assert len(shape['coordinates']) == 1 # no hole - assert val == 127 - try: - shape, val = next(shapes) - except StopIteration: - assert True - else: - assert False - - -def test_shapes_band_shortcut(): - """Test rasterio bands as input to shapes""" - - with rasterio.drivers(): - with rasterio.open('tests/data/shade.tif') as src: - shapes = ftrz.shapes(rasterio.band(src, 1)) - shape, val = next(shapes) - assert shape['type'] == 'Polygon' - assert len(shape['coordinates']) == 1 - assert val == 255 - - -def test_shapes_internal_driver_manager(): - """Make sure this works if driver is managed outside this test""" - - image = numpy.zeros((20, 20), dtype=rasterio.ubyte) - image[5:15, 5:15] = 127 - shapes = ftrz.shapes(image) - shape, val = next(shapes) - assert shape['type'] == 'Polygon' - - -def test_shapes_connectivity(): - """Test connectivity options""" - - image = numpy.zeros((20, 20), dtype=rasterio.ubyte) - image[5:11, 5:11] = 1 - image[11, 11] = 1 - - shapes = ftrz.shapes(image, connectivity=4) - shape, val = next(shapes) - assert len(shape['coordinates'][0]) == 5 - - shapes = ftrz.shapes(image, connectivity=8) - shape, val = next(shapes) - assert len(shape['coordinates'][0]) == 9 - # Note: geometry is not technically valid at this point, it has a self - # intersection at 11,11 - - # Test invalid connectivity - with pytest.raises(ValueError): - shapes = ftrz.shapes(image, connectivity=12) - next(shapes) - - -def test_shapes_dtype(): - """Test image data type handling""" - - rows = cols = 10 - with rasterio.drivers(): - supported_types = ( - ('int16', -32768), - ('int32', -2147483648), - ('uint8', 255), - ('uint16', 65535), - ('float32', 1.434532) - ) - - for dtype, test_value in supported_types: - image = numpy.zeros((rows, cols), dtype=dtype) - image[2:5, 2:5] = test_value - - shapes = ftrz.shapes(image) - shape, value = next(shapes) - if dtype == 'float32': - assert round(value, 6) == round(test_value, 6) - else: - assert value == test_value - - # Unsupported types should all raise exceptions - unsupported_types = ( - ('int8', -127), - ('uint32', 4294967295), - ('int64', 20439845334323), - ('float16', -9343.232), - ('float64', -98332.133422114) - ) - - for dtype, test_value in unsupported_types: - with pytest.raises(ValueError): - image = numpy.zeros((rows, cols), dtype=dtype) - image[2:5, 2:5] = test_value - shapes = ftrz.shapes(image) - next(shapes) - - # Test mask types - image = numpy.zeros((rows, cols), dtype='uint8') - image.fill(255) - supported_mask_types = ( - ('bool', 1), - ('uint8', 255) - ) - for dtype, mask_value in supported_mask_types: - mask = numpy.zeros((rows, cols), dtype=dtype) - mask[2:5, 2:5] = mask_value - shapes = ftrz.shapes(image, mask=mask) - shape, value = next(shapes) - assert value == 255 - - unsupported_mask_types = ( - ('int8', -127), - ('int16', -32768) - ) - for dtype, mask_value in unsupported_mask_types: - with pytest.raises(ValueError): - mask = numpy.zeros((rows, cols), dtype=dtype) - mask[2:5, 2:5] = mask_value - shapes = ftrz.shapes(image, mask=mask) - next(shapes) \ No newline at end of file diff --git a/tests/test_features_sieve.py b/tests/test_features_sieve.py deleted file mode 100644 index e34cf346..00000000 --- a/tests/test_features_sieve.py +++ /dev/null @@ -1,185 +0,0 @@ -import logging -import sys -import numpy -import pytest - -import rasterio -import rasterio.features as ftrz - - -logging.basicConfig(stream=sys.stderr, level=logging.DEBUG) - - -def test_sieve(): - """Test sieving a 10x10 feature from an ndarray.""" - - image = numpy.zeros((20, 20), dtype=rasterio.ubyte) - image[5:15, 5:15] = 1 - - # An attempt to sieve out features smaller than 100 should not change the - # image. - with rasterio.drivers(): - sieved_image = ftrz.sieve(image, 100) - assert numpy.array_equal(sieved_image, image) - - # Setting the size to 100 should leave us an empty, False image. - with rasterio.drivers(): - sieved_image = ftrz.sieve(image, 101) - assert not sieved_image.any() - - # Invalid size value should fail - for invalid_size in (0, 45.1234, image.size + 1): - with pytest.raises(ValueError): - sieved_image = ftrz.sieve(image, invalid_size) - - -def test_sieve_connectivity(): - """Test proper behavior of connectivity""" - - image = numpy.zeros((20, 20), dtype=rasterio.ubyte) - image[5:15:2, 5:15] = 1 - image[6, 4] = 1 - image[8, 15] = 1 - image[10, 4] = 1 - image[12, 15] = 1 - - # Diagonals not connected, all become small features that will be removed - sieved_image = ftrz.sieve(image, 54, connectivity=4) - assert not sieved_image.any() - - # Diagonals connected, everything is retained - sieved_image = ftrz.sieve(image, 54, connectivity=8) - assert numpy.array_equal(sieved_image, image) - - # Invalid connectivity value should fail - with pytest.raises(ValueError): - ftrz.sieve(image, 54, connectivity=12) - - -def test_sieve_output(): - """Test proper behavior of output image, if passed into sieve""" - - with rasterio.drivers(): - shape = (20, 20) - image = numpy.zeros(shape, dtype=rasterio.ubyte) - image[5:15, 5:15] = 1 - - # Output should match returned array - output = numpy.zeros_like(image) - output[1:3, 1:3] = 5 - sieved_image = ftrz.sieve(image, 100, output=output) - assert numpy.array_equal(output, sieved_image) - - # Output of different dtype should fail - output = numpy.zeros(shape, dtype=rasterio.int32) - with pytest.raises(ValueError): - ftrz.sieve(image, 100, output) - - # Output of a different shape should fail - output = numpy.zeros((21, 21), dtype=rasterio.ubyte) - with pytest.raises(ValueError): - ftrz.sieve(image, 100, output) - - -def test_sieve_mask(): - """Test proper behavior of mask image, if passed int sieve""" - - with rasterio.drivers(): - shape = (20, 20) - image = numpy.zeros(shape, dtype=rasterio.ubyte) - image[5:15, 5:15] = 1 - image[1:3, 1:3] = 2 - - # Blank mask has no effect, only areas smaller than size will be - # removed - mask = numpy.ones(shape, dtype=rasterio.bool_) - sieved_image = ftrz.sieve(image, 100, mask=mask) - truth = numpy.zeros_like(image) - truth[5:15, 5:15] = 1 - assert numpy.array_equal(sieved_image, truth) - - # Only areas within the overlap of the mask and values will be kept - mask = numpy.ones(shape, dtype=rasterio.bool_) - mask[7:10, 7:10] = False - sieved_image = ftrz.sieve(image, 100, mask=mask) - truth = numpy.zeros_like(image) - truth[7:10, 7:10] = 1 - assert numpy.array_equal(sieved_image, truth) - - with pytest.raises(ValueError): - mask = numpy.ones((21, 21), dtype=rasterio.bool_) - ftrz.sieve(image, 100, mask=mask) - - -def test_dtypes(): - """Test data type support for sieve""" - - rows = cols = 10 - with rasterio.drivers(): - supported_types = ( - ('int16', -32768), - ('int32', -2147483648), - ('uint8', 255), - ('uint16', 65535) - ) - - for dtype, test_value in supported_types: - image = numpy.zeros((rows, cols), dtype=dtype) - image[2:5, 2:5] = test_value - - # Sieve should return the original image - sieved_image = ftrz.sieve(image, 2) - assert numpy.array_equal(image, sieved_image) - assert numpy.dtype(sieved_image.dtype).name == dtype - - # Sieve should return a blank image - sieved_image = ftrz.sieve(image, 10) - assert numpy.array_equal(numpy.zeros_like(image), sieved_image) - assert numpy.dtype(sieved_image.dtype).name == dtype - - # Unsupported types should all raise exceptions - unsupported_types = ( - ('int8', -127), - ('uint32', 4294967295), - ('int64', 20439845334323), - ('float16', -9343.232), - ('float32', 1.434532), - ('float64', -98332.133422114) - ) - - for dtype, test_value in unsupported_types: - with pytest.raises(ValueError): - image = numpy.zeros((rows, cols), dtype=dtype) - image[2:5, 2:5] = test_value - ftrz.sieve(image, 2) - - # Test mask types - image = numpy.zeros((rows, cols), dtype='uint8') - image.fill(255) - - supported_mask_types = ( - ('bool', 1), - ('uint8', 255) - ) - for dtype, mask_value in supported_mask_types: - mask = numpy.zeros((rows, cols), dtype=dtype) - mask[2:5, 2:5] = mask_value - sieved_image = ftrz.sieve(image, 2, mask=mask) - assert numpy.array_equal(image, sieved_image) - - unsupported_mask_types = ( - ('int8', -127), - ('int16', -32768) - ) - for dtype, mask_value in unsupported_mask_types: - with pytest.raises(ValueError): - mask = numpy.zeros((rows, cols), dtype=dtype) - mask[2:5, 2:5] = mask_value - ftrz.sieve(image, 2, mask=mask) - - -def test_sieve_shade(): - with rasterio.drivers(): - with rasterio.open('tests/data/shade.tif') as src: - sieved_image = ftrz.sieve(rasterio.band(src, 1), 42) - assert sieved_image.shape == (1024, 1024) diff --git a/tests/test_rio_features.py b/tests/test_rio_features.py index a1a4b153..b1f20a0c 100644 --- a/tests/test_rio_features.py +++ b/tests/test_rio_features.py @@ -3,559 +3,738 @@ import os import re import sys import numpy +import json +from affine import Affine import rasterio from rasterio.rio import features +DEFAULT_SHAPE = (10, 10) + logging.basicConfig(stream=sys.stderr, level=logging.DEBUG) -TEST_FEATURES = """{ - "geometry": { - "coordinates": [ - [ - [-110, 40], - [-100, 40], - [-100, 45], - [-105, 45], - [-110, 40] - ] - ], - "type": "Polygon" - }, - "properties": { - "val": 15 - }, - "type": "Feature" -}""" -TEST_MERC_FEATURES = """{ - "geometry": { - "coordinates": [ - [ - [-11858134, 4808920], - [-11868134, 4804143], - [-11853357, 4804143], - [-11840000, 4812000], - [-11858134, 4808920] - ] - ], - "type": "Polygon" - }, - "properties": { - "val": 10 - }, - "type": "Feature" -}""" +def test_mask(runner, tmpdir, basic_feature, basic_image_2x2, + pixelated_image_file): -# > rio shapes tests/data/shade.tif --mask --sampling 500 --projected --precision 0 -TEST_MERC_FEATURECOLLECTION = """{ - "bbox": [-11858135.0, 4803914.0, -11848351.0, 4813698.0], - "features": [{ - "bbox": [-11853357.504145855, 4808920.97837715, - -11848580.189878704, 4813698.2926443005], - "geometry": { - "coordinates": [ - [ - [-11853357.504145855, 4813698.2926443005], - [-11853357.504145855, 4808920.97837715], - [-11848580.189878704, 4808920.97837715], - [-11848580.189878704, 4813698.2926443005], - [-11853357.504145855, 4813698.2926443005] - ] - ], - "type": "Polygon" - }, - "properties": { - "val": 2 - }, - "type": "Feature" - }, { - "bbox": [-11858134.818413004, 4804143.66411, - -11853357.504145855, 4808920.97837715], - "geometry": { - "coordinates": [ - [ - [-11858134.818413004, 4808920.97837715], - [-11858134.818413004, 4804143.66411], - [-11853357.504145855, 4804143.66411], - [-11853357.504145855, 4808920.97837715], - [-11858134.818413004, 4808920.97837715] - ] - ], - "type": "Polygon" - }, - "properties": { - "val": 3 - }, - "type": "Feature" - }], - "type": "FeatureCollection" -}""" - - -def test_mask(runner, tmpdir): output = str(tmpdir.join('test.tif')) - with rasterio.open('tests/data/shade.tif') as src: - src_data = src.read(1, masked=True) - - result = runner.invoke( - features.mask, - ['tests/data/shade.tif', output, '--geojson-mask', '-'], - input=TEST_MERC_FEATURES - ) - assert result.exit_code == 0 - assert os.path.exists(output) - - masked_count = 0 - with rasterio.open(output) as out: - assert out.count == src.count - assert out.shape == src.shape - out_data = out.read(1, masked=True) - - # Make sure that pixels with 0 value were converted to mask - masked_count = (src_data == 0).sum() - (out_data == 0).sum() - assert masked_count == 79743 - assert out_data.mask.sum() - src_data.mask.sum() == masked_count - - # Test using --all-touched option - result = runner.invoke( - features.mask, - [ - 'tests/data/shade.tif', output, - '--all', - '--geojson-mask', '-' - ], - input=TEST_MERC_FEATURES - ) - assert result.exit_code == 0 - with rasterio.open(output) as out: - out_data = out.read(1, masked=True) - - # Make sure that more pixels with 0 value were converted to mask - masked_count2 = (src_data == 0).sum() - (out_data == 0).sum() - assert masked_count2 > 0 and masked_count > masked_count2 - assert out_data.mask.sum() - src_data.mask.sum() == masked_count2 - - # Test using --invert option - result = runner.invoke( - features.mask, - [ - 'tests/data/shade.tif', output, - '--invert', - '--geojson-mask', '-' - ], - input=TEST_MERC_FEATURES - ) - assert result.exit_code == 0 - with rasterio.open(output) as out: - out_data = out.read(1, masked=True) - # Areas that were masked when not inverted should now be 0 - assert (out_data == 0).sum() == masked_count - - # Test with feature collection result = runner.invoke( features.mask, - ['tests/data/shade.tif', output, '--geojson-mask', '-'], - input=TEST_MERC_FEATURECOLLECTION + [pixelated_image_file, output, '--geojson-mask', '-'], + input=json.dumps(basic_feature) + ) + + assert result.exit_code == 0 + assert os.path.exists(output) + + with rasterio.open(output) as out: + assert numpy.array_equal( + basic_image_2x2, + out.read(1, masked=True).filled(0) + ) + + +def test_mask_all_touched(runner, tmpdir, basic_feature, basic_image, + pixelated_image_file): + + output = str(tmpdir.join('test.tif')) + + result = runner.invoke( + features.mask, + [pixelated_image_file, output, '--all', '--geojson-mask', '-'], + input=json.dumps(basic_feature) ) assert result.exit_code == 0 + assert os.path.exists(output) + + with rasterio.open(output) as out: + assert numpy.array_equal( + basic_image, + out.read(1, masked=True).filled(0) + ) + + +def test_mask_invert(runner, tmpdir, basic_feature, pixelated_image, + pixelated_image_file): + + truth = pixelated_image + truth[2:4, 2:4] = 0 + + output = str(tmpdir.join('test.tif')) - # Missing GeoJSON should make copy of input to output - output2 = str(tmpdir.join('test2.tif')) result = runner.invoke( features.mask, - ['tests/data/shade.tif', output2] + [pixelated_image_file, output, '--invert', '--geojson-mask', '-'], + input=json.dumps(basic_feature) ) assert result.exit_code == 0 - assert os.path.exists(output2) - with rasterio.open('tests/data/shade.tif') as src: - with rasterio.open(output2) as out: - src_data = src.read(1, masked=True) - out_data = out.read(1, masked=True) - assert numpy.array_equal(src_data, out_data) + assert os.path.exists(output) + + with rasterio.open(output) as out: + assert numpy.array_equal( + truth, + out.read(1, masked=True).filled(0) + ) + + +def test_mask_featurecollection(runner, tmpdir, basic_featurecollection, + basic_image_2x2, pixelated_image_file): + + output = str(tmpdir.join('test.tif')) - # Invalid JSON should fail result = runner.invoke( features.mask, - ['tests/data/shade.tif', output, '--geojson-mask', '-'], + [pixelated_image_file, output, '--geojson-mask', '-'], + input=json.dumps(basic_featurecollection) + ) + assert result.exit_code == 0 + assert os.path.exists(output) + + with rasterio.open(output) as out: + assert numpy.array_equal( + basic_image_2x2, + out.read(1, masked=True).filled(0) + ) + + +def test_mask_out_of_bounds(runner, tmpdir, basic_feature, + pixelated_image_file): + """ + A GeoJSON mask that is outside bounds of raster should result in a + blank image. + """ + + coords = numpy.array(basic_feature['geometry']['coordinates']) - 10 + basic_feature['geometry']['coordinates'] = coords.tolist() + + output = str(tmpdir.join('test.tif')) + + result = runner.invoke( + features.mask, + [pixelated_image_file, output, '--geojson-mask', '-'], + input=json.dumps(basic_feature) + ) + assert result.exit_code == 0 + assert 'outside bounds' in result.output + assert os.path.exists(output) + + with rasterio.open(output) as out: + assert not numpy.any(out.read(1, masked=True).filled(0)) + + +def test_mask_no_geojson(runner, tmpdir, pixelated_image, pixelated_image_file): + """ Mask without geojson input should simply return same raster as input """ + + output = str(tmpdir.join('test.tif')) + + result = runner.invoke(features.mask, [pixelated_image_file, output]) + assert result.exit_code == 0 + assert os.path.exists(output) + + with rasterio.open(output) as out: + assert numpy.array_equal( + pixelated_image, + out.read(1, masked=True).filled(0) + ) + + +def test_mask_invalid_geojson(runner, tmpdir, pixelated_image_file): + """ Invalid GeoJSON should fail """ + + output = str(tmpdir.join('test.tif')) + + # Using invalid JSON + result = runner.invoke( + features.mask, + [pixelated_image_file, output, '--geojson-mask', '-'], input='{bogus: value}' ) assert result.exit_code == 2 assert 'GeoJSON could not be read' in result.output + # Using invalid GeoJSON result = runner.invoke( features.mask, - ['tests/data/shade.tif', output, '--geojson-mask', '-'], + [pixelated_image_file, output, '--geojson-mask', '-'], input='{"bogus": "value"}' ) assert result.exit_code == 2 assert 'Invalid GeoJSON' in result.output -def test_mask_crop(runner, tmpdir): +def test_mask_crop(runner, tmpdir, basic_feature, pixelated_image): + """ + In order to test --crop option, we need to use a transform more similar to + a normal raster, with a negative y pixel size. + """ + + image = pixelated_image + outfilename = str(tmpdir.join('pixelated_image.tif')) + kwargs = { + "crs": {'init': 'epsg:4326'}, + "transform": Affine(1, 0, 0, 0, -1, 0), + "count": 1, + "dtype": rasterio.uint8, + "driver": "GTiff", + "width": image.shape[1], + "height": image.shape[0], + "nodata": 255 + } + with rasterio.drivers(): + with rasterio.open(outfilename, 'w', **kwargs) as out: + out.write_band(1, image) + + output = str(tmpdir.join('test.tif')) - with rasterio.open('tests/data/shade.tif') as src: + truth = numpy.zeros((4, 3)) + truth[1:3, 0:2] = 1 - result = runner.invoke( - features.mask, - [ - 'tests/data/shade.tif', output, - '--crop', - '--geojson-mask', '-' - ], - input=TEST_MERC_FEATURES + result = runner.invoke( + features.mask, + [outfilename, output, '--crop', '--geojson-mask', '-'], + input=json.dumps(basic_feature) + ) + + assert result.exit_code == 0 + assert os.path.exists(output) + with rasterio.open(output) as out: + assert numpy.array_equal( + truth, + out.read(1, masked=True).filled(0) ) - assert result.exit_code == 0 - assert os.path.exists(output) - with rasterio.open(output) as out: - assert out.shape[1] == src.shape[1] - assert out.shape[0] < src.shape[0] - assert out.shape[0] == 824 - - # Adding invert option after crop should be ignored - result = runner.invoke( - features.mask, - [ - 'tests/data/shade.tif', output, - '--crop', - '--invert', - '--geojson-mask', '-' - ], - input=TEST_MERC_FEATURES - ) - assert result.exit_code == 0 - assert 'Invert option ignored' in result.output -def test_mask_out_of_bounds(runner, tmpdir): +def test_mask_crop_inverted_y(runner, tmpdir, basic_feature, pixelated_image_file): + """ + --crop option should also work if raster has a positive y pixel size + (e.g., Affine.identity() ). + """ + output = str(tmpdir.join('test.tif')) - # Crop with out of bounds raster should - result = runner.invoke( - features.mask, - ['tests/data/shade.tif', output, '--geojson-mask', '-'], - input=TEST_FEATURES - ) - assert result.exit_code == 0 - assert 'outside bounds' in result.output - # Crop with out of bounds raster should fail + truth = numpy.zeros((4, 3)) + truth[1:3, 0:2] = 1 + result = runner.invoke( features.mask, - [ - 'tests/data/shade.tif', output, - '--crop', - '--geojson-mask', '-' - ], - input=TEST_FEATURES + [pixelated_image_file, output, '--crop', '--geojson-mask', '-'], + input=json.dumps(basic_feature) + ) + + assert result.exit_code == 0 + assert os.path.exists(output) + with rasterio.open(output) as out: + assert numpy.array_equal( + truth, + out.read(1, masked=True).filled(0) + ) + + +def test_mask_crop_out_of_bounds(runner, tmpdir, basic_feature, + pixelated_image_file): + """ + A GeoJSON mask that is outside bounds of raster should fail with + --crop option. + """ + + coords = numpy.array(basic_feature['geometry']['coordinates']) - 10 + basic_feature['geometry']['coordinates'] = coords.tolist() + + output = str(tmpdir.join('test.tif')) + + result = runner.invoke( + features.mask, + [pixelated_image_file, output, '--crop', '--geojson-mask', '-'], + input=json.dumps(basic_feature) ) assert result.exit_code == 2 assert 'not allowed' in result.output -def test_shapes(runner): - result = runner.invoke(features.shapes, ['tests/data/shade.tif']) +def test_mask_crop_and_invert(runner, tmpdir, basic_feature, pixelated_image, + pixelated_image_file): + """ Adding crop and invert options should ignore invert option """ + + output = str(tmpdir.join('test.tif')) + + result = runner.invoke( + features.mask, + [ + pixelated_image_file, output, + '--crop', + '--invert', + '--geojson-mask', '-' + ], + input=json.dumps(basic_feature) + ) + assert result.exit_code == 0 + assert 'Invert option ignored' in result.output + + +def test_shapes(runner, pixelated_image_file): + result = runner.invoke(features.shapes, [pixelated_image_file]) + assert result.exit_code == 0 assert result.output.count('"FeatureCollection"') == 1 - assert result.output.count('"Feature"') == 232 + assert result.output.count('"Feature"') == 4 + assert numpy.allclose( + json.loads(result.output)['features'][0]['geometry']['coordinates'], + [[[0, 0], [0, 1], [1, 1], [1, 0], [0, 0]]] + ) + + +def test_shapes_invalid_bidx(runner, pixelated_image_file): + result = runner.invoke(features.shapes, [pixelated_image_file, '--bidx', 4]) - # Invalid band index should fail - result = runner.invoke( - features.shapes, ['tests/data/shade.tif', '--bidx', '4']) assert result.exit_code == 1 + # Underlying exception message trapped by shapes -def test_shapes_sequence(runner): - result = runner.invoke(features.shapes, ['tests/data/shade.tif', '--sequence']) +def test_shapes_sequence(runner, pixelated_image_file): + """ + --sequence option should produce 4 features in series rather than + inside a feature collection. + """ + + result = runner.invoke(features.shapes, [pixelated_image_file, '--sequence']) + assert result.exit_code == 0 assert result.output.count('"FeatureCollection"') == 0 - assert result.output.count('"Feature"') == 232 + assert result.output.count('"Feature"') == 4 + assert result.output.count('\n') == 4 -def test_shapes_sequence_rs(runner): +def test_shapes_sequence_rs(runner, pixelated_image_file): + """ --rs option should use the feature separator character. """ + result = runner.invoke( - features.shapes, [ - 'tests/data/shade.tif', - '--sequence', - '--rs']) + features.shapes, [pixelated_image_file, '--sequence', '--rs'] + ) + assert result.exit_code == 0 assert result.output.count('"FeatureCollection"') == 0 - assert result.output.count('"Feature"') == 232 - assert result.output.count(u'\u001e') == 232 + assert result.output.count('"Feature"') == 4 + assert result.output.count(u'\u001e') == 4 -def test_shapes_with_nodata(runner): - result = runner.invoke(features.shapes, ['tests/data/shade.tif', '--with-nodata']) + +def test_shapes_with_nodata(runner, pixelated_image, pixelated_image_file): + """ + An area of nodata should also be represented with a shape when using + --with-nodata option + """ + + pixelated_image[0:2, 8:10] = 255 + + with rasterio.open(pixelated_image_file, 'r+') as out: + out.write_band(1, pixelated_image) + + result = runner.invoke( + features.shapes, [pixelated_image_file, '--with-nodata'] + ) assert result.exit_code == 0 assert result.output.count('"FeatureCollection"') == 1 - assert result.output.count('"Feature"') == 288 + assert result.output.count('"Feature"') == 5 -def test_shapes_indent(runner): - result = runner.invoke(features.shapes, ['tests/data/shade.tif', '--indent', '2']) +def test_shapes_indent(runner, pixelated_image_file): + """ + --indent option should produce lots of newlines and contiguous spaces + """ + + result = runner.invoke( + features.shapes, [pixelated_image_file, '--indent', 2] + ) + assert result.exit_code == 0 assert result.output.count('"FeatureCollection"') == 1 - assert result.output.count('\n') == 70371 + assert result.output.count('"Feature"') == 4 + assert result.output.count('\n') == 231 + assert result.output.count(' ') == 180 -def test_shapes_compact(runner): - result = runner.invoke(features.shapes, ['tests/data/shade.tif', '--compact']) +def test_shapes_compact(runner, pixelated_image_file): + result = runner.invoke(features.shapes, [pixelated_image_file, '--compact']) + assert result.exit_code == 0 assert result.output.count('"FeatureCollection"') == 1 + assert result.output.count('"Feature"') == 4 assert result.output.count(', ') == 0 assert result.output.count(': ') == 0 -def test_shapes_sampling(runner): +def test_shapes_sampling(runner, pixelated_image_file): + """ --sampling option should remove the single pixel features """ result = runner.invoke( - features.shapes, ['tests/data/shade.tif', '--sampling', '11']) + features.shapes, [pixelated_image_file, '--sampling', 2] + ) + assert result.exit_code == 0 assert result.output.count('"FeatureCollection"') == 1 - assert result.output.count('"Feature"') == 117 + assert result.output.count('"Feature"') == 2 -def test_shapes_precision(runner): +def test_shapes_precision(runner, pixelated_image_file): + """ Output numbers should have no more than 1 decimal place """ + result = runner.invoke( - features.shapes, ['tests/data/shade.tif', '--precision', '1']) + features.shapes, [pixelated_image_file, '--precision', 1] + ) + assert result.exit_code == 0 assert result.output.count('"FeatureCollection"') == 1 - # Find no numbers with 2+ decimal places. + assert result.output.count('"Feature"') == 4 assert re.search(r'\d*\.\d{2,}', result.output) is None -def test_shapes_mask(runner): - result = runner.invoke(features.shapes, ['tests/data/RGB.byte.tif', '--mask']) - assert result.exit_code == 0 - assert result.output.count('"FeatureCollection"') == 1 - assert result.output.count('"Feature"') == 7 +def test_shapes_mask(runner, pixelated_image, pixelated_image_file): + """ --mask should extract the nodata area of the image """ + pixelated_image[0:5, 0:10] = 255 + pixelated_image[0:10, 0:3] = 255 + pixelated_image[8:10, 8:10] = 255 + + with rasterio.open(pixelated_image_file, 'r+') as out: + out.write_band(1, pixelated_image) + + result = runner.invoke(features.shapes, [pixelated_image_file, '--mask']) + + print(result.output) + print(result.exception) -def test_shapes_mask_decimated(runner): - result = runner.invoke( - features.shapes, - ['tests/data/RGB.byte.tif', '--mask', '--sampling', '10']) assert result.exit_code == 0 assert result.output.count('"FeatureCollection"') == 1 assert result.output.count('"Feature"') == 1 + assert numpy.allclose( + json.loads(result.output)['features'][0]['geometry']['coordinates'], + [[[3, 5], [3, 10], [8, 10], [8, 8], [9, 8], [10, 8], [10, 5], [3, 5]]] + ) + + +def test_shapes_mask_sampling(runner, pixelated_image, pixelated_image_file): + """ + using --sampling with the mask should snap coordinates to the nearest + factor of 5 + """ + pixelated_image[0:5, 0:10] = 255 + pixelated_image[0:10, 0:3] = 255 + pixelated_image[8:10, 8:10] = 255 + + with rasterio.open(pixelated_image_file, 'r+') as out: + out.write_band(1, pixelated_image) + + result = runner.invoke( + features.shapes, [pixelated_image_file, '--mask', '--sampling', 5] + ) -def test_shapes_band1_as_mask(runner): - result = runner.invoke(features.shapes, - ['tests/data/RGB.byte.tif', '--band', '--bidx', '1', '--as-mask']) assert result.exit_code == 0 assert result.output.count('"FeatureCollection"') == 1 - assert result.output.count('"Feature"') == 9 + assert result.output.count('"Feature"') == 1 + + assert numpy.allclose( + json.loads(result.output)['features'][0]['geometry']['coordinates'], + [[[5, 5], [5, 10], [10, 10], [10, 5], [5, 5]]] + ) -def test_rasterize_err(tmpdir, runner): +def test_shapes_band1_as_mask(runner, pixelated_image, pixelated_image_file): + """ + When using --as-mask option, pixel value should not matter, only depends + on pixels being contiguous. + """ + + pixelated_image[2:3, 2:3] = 4 + + with rasterio.open(pixelated_image_file, 'r+') as out: + out.write_band(1, pixelated_image) + + result = runner.invoke( + features.shapes, + [pixelated_image_file, '--band', '--bidx', '1', '--as-mask'] + ) + + assert result.exit_code == 0 + assert result.output.count('"FeatureCollection"') == 1 + assert result.output.count('"Feature"') == 3 + assert numpy.allclose( + json.loads(result.output)['features'][1]['geometry']['coordinates'], + [[[2, 2], [2, 5], [5, 5], [5, 2], [2, 2]]] + ) + + +def test_rasterize(tmpdir, runner, basic_feature): output = str(tmpdir.join('test.tif')) - # Test invalid stdin - result = runner.invoke(features.rasterize, [output], input='BOGUS') - assert result.exit_code == -1 + result = runner.invoke( + features.rasterize, + [output, '--dimensions', DEFAULT_SHAPE[0], DEFAULT_SHAPE[1]], + input=json.dumps(basic_feature) + ) - # Test invalid GeoJSON - result = runner.invoke(features.rasterize, [output], - input='{"foo": "bar"}') - assert result.exit_code == 2 - - # Test invalid res - result = runner.invoke(features.rasterize, [output], input=TEST_FEATURES) - assert result.exit_code == 2 - - # Test invalid CRS for bounds - result = runner.invoke(features.rasterize, [output, '--res', 1], - input=TEST_MERC_FEATURECOLLECTION) - assert result.exit_code == 2 - - # Test invalid CRS value - result = runner.invoke(features.rasterize, [output, - '--res', 1, - '--src-crs', 'BOGUS'], - input=TEST_MERC_FEATURECOLLECTION) - assert result.exit_code == 2 + assert result.exit_code == 0 + assert os.path.exists(output) + with rasterio.open(output) as out: + assert numpy.allclose(out.bounds, (2, 2, 4.25, 4.25)) + data = out.read(1, masked=False) + assert data.shape == DEFAULT_SHAPE + assert numpy.all(data) -def test_rasterize(tmpdir, runner): - # Test dimensions +def test_rasterize_bounds(tmpdir, runner, basic_feature, basic_image_2x2): output = str(tmpdir.join('test.tif')) - result = runner.invoke(features.rasterize, - [output, '--dimensions', 20, 10], - input=TEST_FEATURES) + result = runner.invoke( + features.rasterize, + [ + output, + '--dimensions', DEFAULT_SHAPE[0], DEFAULT_SHAPE[1], + '--bounds', 0, 10, 10, 0 + ], + input=json.dumps(basic_feature) + ) + assert result.exit_code == 0 assert os.path.exists(output) with rasterio.open(output) as out: - assert out.count == 1 - assert out.meta['width'] == 20 - assert out.meta['height'] == 10 + assert numpy.allclose(out.bounds, (0, 10, 10, 0)) data = out.read(1, masked=False) - assert (data == 0).sum() == 55 - assert (data == 1).sum() == 145 + assert numpy.array_equal(basic_image_2x2, data) + assert data.shape == DEFAULT_SHAPE + + +def test_rasterize_resolution(tmpdir, runner, basic_feature): + output = str(tmpdir.join('test.tif')) + result = runner.invoke( + features.rasterize, + [output, '--res', 0.15], + input=json.dumps(basic_feature) + ) - # Test dimensions and bounds - output = str(tmpdir.join('test2.tif')) - result = runner.invoke(features.rasterize, - [output, - '--dimensions', 40, 20, - '--bounds', -120, 30, -90, 50 - ], input=TEST_FEATURES) assert result.exit_code == 0 assert os.path.exists(output) with rasterio.open(output) as out: - assert out.count == 1 - assert out.meta['width'] == 40 - assert out.meta['height'] == 20 + assert numpy.allclose(out.bounds, (2, 2, 4.25, 4.25)) data = out.read(1, masked=False) - assert (data == 0).sum() == 748 - assert (data == 1).sum() == 52 + assert data.shape == (15, 15) + assert numpy.all(data) + + +def test_rasterize_src_crs(tmpdir, runner, basic_feature): + output = str(tmpdir.join('test.tif')) + result = runner.invoke( + features.rasterize, + [ + output, + '--dimensions', DEFAULT_SHAPE[0], DEFAULT_SHAPE[1], + '--src-crs', 'EPSG:3857' + ], + input=json.dumps(basic_feature) + ) - # Test resolution - output = str(tmpdir.join('test3.tif')) - result = runner.invoke(features.rasterize, - [output, '--res', 0.5], input=TEST_FEATURES) assert result.exit_code == 0 assert os.path.exists(output) - with rasterio.open(output) as out: - assert out.count == 1 - assert out.meta['width'] == 20 - assert out.meta['height'] == 10 - data = out.read(1, masked=False) - assert (data == 0).sum() == 55 - assert (data == 1).sum() == 145 - - # Test that src-crs is written into new output - output = str(tmpdir.join('test4.tif')) - result = runner.invoke(features.rasterize, - [output, - '--dimensions', 20, 10, - '--src-crs', 'EPSG:3857' - ], - input=TEST_MERC_FEATURECOLLECTION) - assert result.exit_code == 0 with rasterio.open(output) as out: assert out.crs['init'].lower() == 'epsg:3857' -def test_rasterize_existing_output(tmpdir, runner): +def test_rasterize_mismatched_src_crs(tmpdir, runner, basic_feature): + """ + A --src-crs that is geographic with coordinates that are outside + world bounds should fail. + """ + + coords = numpy.array(basic_feature['geometry']['coordinates']) * 100000 + basic_feature['geometry']['coordinates'] = coords.tolist() + output = str(tmpdir.join('test.tif')) - result = runner.invoke(features.rasterize, - [output, '--res', 0.5], input=TEST_FEATURES) - assert result.exit_code == 0 - assert os.path.exists(output) + result = runner.invoke( + features.rasterize, + [ + output, + '--dimensions', DEFAULT_SHAPE[0], DEFAULT_SHAPE[1], + '--src-crs', 'EPSG:4326' + ], + input=json.dumps(basic_feature) + ) - geojson = """{ - "geometry": { - "coordinates": [ - [ - [-102, 40], - [-98, 40], - [-98, 45], - [-100, 45], - [-102, 40] - ] - ], - "type": "Polygon" - }, - "type": "Feature" - }""" - - result = runner.invoke(features.rasterize, [output, '--default-value', 2], - input=geojson) - - with rasterio.open(output) as out: - assert out.count == 1 - data = out.read(1, masked=False) - assert (data == 0).sum() == 55 - assert (data == 1).sum() == 125 - assert (data == 2).sum() == 20 - - # Confirm that a different src-crs is rejected, even if a geographic crs - result = runner.invoke(features.rasterize, - [output, - '--res', 0.5, - '--src-crs', 'EPSG:4269' - ], input=TEST_FEATURES) assert result.exit_code == 2 + assert 'Bounds are beyond the valid extent for EPSG:4326' in result.output -def test_rasterize_like(tmpdir, runner): +def test_rasterize_invalid_src_crs(tmpdir, runner, basic_feature): output = str(tmpdir.join('test.tif')) - result = runner.invoke(features.rasterize, - [output, '--like', 'tests/data/shade.tif'], - input=TEST_MERC_FEATURECOLLECTION) - assert result.exit_code == 0 - assert os.path.exists(output) - with rasterio.open(output) as out: - assert out.count == 1 - data = out.read(1, masked=False) - assert (data == 0).sum() == 548576 - assert (data == 1).sum() == 500000 + result = runner.invoke( + features.rasterize, + [ + output, + '--dimensions', DEFAULT_SHAPE[0], DEFAULT_SHAPE[1], + '--src-crs', 'foo:bar' + ], + input=json.dumps(basic_feature) + ) - # Test invalid like raster - output = str(tmpdir.join('test2.tif')) - result = runner.invoke(features.rasterize, - [output, '--like', str(tmpdir.join('foo.tif'))], input=TEST_FEATURES) assert result.exit_code == 2 + assert 'invalid CRS. Must be an EPSG code.' in result.output + + +def test_rasterize_existing_output(tmpdir, runner, basic_feature): + """ + Create a rasterized output, then rasterize additional pixels into it. + The final result should include rasterized pixels from both features. + """ + + truth = numpy.zeros(DEFAULT_SHAPE) + truth[2:4, 2:4] = 1 + truth[4:6, 4:6] = 1 + + output = str(tmpdir.join('test.tif')) + result = runner.invoke( + features.rasterize, + [ + output, + '--dimensions', DEFAULT_SHAPE[0], DEFAULT_SHAPE[1], + '--bounds', 0, 10, 10, 0 + ], + input=json.dumps(basic_feature) + ) + + assert result.exit_code == 0 + assert os.path.exists(output) + + coords = numpy.array(basic_feature['geometry']['coordinates']) + 2 + basic_feature['geometry']['coordinates'] = coords.tolist() + + result = runner.invoke( + features.rasterize, + [output, '--dimensions', DEFAULT_SHAPE[0], DEFAULT_SHAPE[1]], + input=json.dumps(basic_feature) + ) + + assert result.exit_code == 0 + + with rasterio.open(output) as out: + assert numpy.array_equal(truth, out.read(1, masked=False)) + + +def test_rasterize_like_raster(tmpdir, runner, basic_feature, basic_image_2x2, + pixelated_image_file): + + output = str(tmpdir.join('test.tif')) + + result = runner.invoke( + features.rasterize, + [output, '--like', pixelated_image_file], + input=json.dumps(basic_feature) + ) + + assert result.exit_code == 0 + assert os.path.exists(output) + with rasterio.open(output) as out: + assert numpy.array_equal(basic_image_2x2, out.read(1, masked=False)) + + with rasterio.open(pixelated_image_file) as src: + assert out.crs == src.crs + assert out.bounds == src.bounds + assert src.affine == src.affine + + +def test_rasterize_invalid_like_raster(tmpdir, runner, basic_feature): + output = str(tmpdir.join('test.tif')) + result = runner.invoke( + features.rasterize, + [output, '--like', str(tmpdir.join('foo.tif'))], + input=json.dumps(basic_feature) + ) - # Test that src-crs different than --like raster crs breaks - output = str(tmpdir.join('test3.tif')) - result = runner.invoke(features.rasterize, - [output, - '--like', 'tests/data/shade.tif', - '--src-crs', 'EPSG:4326'], - input=TEST_FEATURES) assert result.exit_code == 2 + assert 'Invalid value for "--like":' in result.output -def test_rasterize_property_value(tmpdir, runner): - # Test feature collection property values +def test_rasterize_like_raster_src_crs_mismatch(tmpdir, runner, basic_feature, + pixelated_image_file): output = str(tmpdir.join('test.tif')) - result = runner.invoke(features.rasterize, - [output, - '--res', 1000, - '--property', 'val', - '--src-crs', 'EPSG:3857' - ], - input=TEST_MERC_FEATURECOLLECTION) + result = runner.invoke( + features.rasterize, + [output, '--like', pixelated_image_file, '--src-crs', 'EPSG:3857'], + input=json.dumps(basic_feature) + ) + + assert result.exit_code == 2 + assert 'GeoJSON does not match crs of --like raster' in result.output + + +def test_rasterize_property_value(tmpdir, runner, basic_feature): + output = str(tmpdir.join('test.tif')) + result = runner.invoke( + features.rasterize, + [ + output, + '--dimensions', DEFAULT_SHAPE[0], DEFAULT_SHAPE[1], + '--property', 'val' + ], + input=json.dumps(basic_feature) + ) + assert result.exit_code == 0 assert os.path.exists(output) with rasterio.open(output) as out: - assert out.count == 1 + assert numpy.allclose(out.bounds, (2, 2, 4.25, 4.25)) data = out.read(1, masked=False) - assert (data == 0).sum() == 50 - assert (data == 2).sum() == 25 - assert (data == 3).sum() == 25 - - # Test feature property values - output = str(tmpdir.join('test2.tif')) - result = runner.invoke(features.rasterize, - [output, '--res', 0.5, '--property', 'val'], - input=TEST_FEATURES) - assert result.exit_code == 0 - assert os.path.exists(output) - with rasterio.open(output) as out: - assert out.count == 1 - data = out.read(1, masked=False) - assert (data == 0).sum() == 55 - assert (data == 15).sum() == 145 + assert data.shape == DEFAULT_SHAPE + assert numpy.all(data == basic_feature['properties']['val']) -def test_rasterize_out_of_bounds(tmpdir, runner): +def test_rasterize_like_raster_outside_bounds(tmpdir, runner, basic_feature, + pixelated_image_file): + """ + Rasterizing a feature outside bounds of --like raster should result + in a blank image + """ + + coords = numpy.array(basic_feature['geometry']['coordinates']) + 100 + basic_feature['geometry']['coordinates'] = coords.tolist() + output = str(tmpdir.join('test.tif')) + result = runner.invoke( + features.rasterize, + [output, '--like', pixelated_image_file], + input=json.dumps(basic_feature) + ) - # Test out of bounds of --like raster - result = runner.invoke(features.rasterize, - [output, '--like', 'tests/data/shade.tif'], - input=TEST_FEATURES) assert result.exit_code == 0 assert 'outside bounds' in result.output assert os.path.exists(output) with rasterio.open(output) as out: - data = out.read_band(1, masked=False) - assert data.sum() == 0 + assert not numpy.any(out.read_band(1, masked=False)) - # Confirm that this does not fail when out of bounds for existing raster - result = runner.invoke(features.rasterize, [output], input=TEST_FEATURES) - assert result.exit_code == 0 - assert 'outside bounds' in result.output + +def test_rasterize_invalid_stdin(tmpdir, runner): + """ Invalid value for stdin should fail with exception """ + + output = str(tmpdir.join('test.tif')) + result = runner.invoke(features.rasterize, [output], input='BOGUS') + + assert result.exit_code == -1 + + +def test_rasterize_invalid_geojson(tmpdir, runner): + """ Invalid GeoJSON should fail with error """ + output = str(tmpdir.join('test.tif')) + result = runner.invoke(features.rasterize, [output], input='{"A": "B"}') + + assert result.exit_code == 2 + assert 'Invalid GeoJSON' in result.output + + +def test_rasterize_missing_parameters(tmpdir, runner, basic_feature): + """ At least --res or --dimensions are required """ + + output = str(tmpdir.join('test.tif')) + result = runner.invoke( + features.rasterize, [output], input=json.dumps(basic_feature) + ) + + assert result.exit_code == 2 + assert 'pixel dimensions are required' in result.output