Merge tag '1.2.4'

This commit is contained in:
Sean Gillies 2021-05-31 16:26:51 -06:00
commit 2de62d2ead
26 changed files with 396 additions and 210 deletions

View File

@ -48,11 +48,11 @@ jobs:
- python: "3.8"
env:
GDALVERSION="master"
PROJVERSION="7.2.1"
PROJVERSION="8.0.1"
allow_failures:
- env:
GDALVERSION="master"
PROJVERSION="7.2.1"
PROJVERSION="8.0.1"
addons:
apt:

View File

@ -10,6 +10,28 @@ Changes
environment variable.
- Add version to CRS.to_wkt() for WKT2 support (#2122).
1.2.4 (2021-05-31)
------------------
- Eliminate unneeded marker for CLI nodata options to be ignored. We will stick
with None (#2191).
- Guard against use of gauss resampling when creating a WarpedVRT (#2190).
- Prevent segfaults when buffer and band indexes are mismatched in
io_multi_band and io_multi_mask (#2189).
- Change comparisons of nodata to IgnoreOption to accomodate changes in click
8.0.
- Use Window constructor instead of from_slices in windows.union to allow a
proper union to be formed from windows extending outside a dataset (#2186).
- Read GDAL CInt16 data as np.complex64 and allow saving complex data to CInt16
(#2185).
- Skip merge sources which do not overlap the destination.
- Allow unsafe casting in the merge tool, restoring behavior of version 1.2.0
(#2179).
- Remove a workaround for an old Python 3.4 bug from the BoundingBox
implementation (#2172).
- Add setuptools as an explicit installation requirement.
>>>>>>> 1.2.4
1.2.3 (2021-04-26)
------------------

View File

@ -171,6 +171,7 @@ cdef _band_dtype(GDALRasterBandH band):
else:
return 'uint8'
return dtypes.dtype_fwd[gdal_dtype]
@ -549,7 +550,9 @@ cdef class DatasetBase:
# to check that the return value is within the range of the
# data type. If so, the band has a nodata value. If not,
# there's no nodata value.
if (success == 0 or
if dtype not in dtypes.dtype_ranges:
pass
elif (success == 0 or
val < dtypes.dtype_ranges[dtype][0] or
val > dtypes.dtype_ranges[dtype][1]):
val = None

View File

@ -9,6 +9,7 @@ import logging
import numpy as np
from rasterio import dtypes
from rasterio.dtypes import _getnpdtype
from rasterio.enums import MergeAlg
from rasterio._err cimport exc_wrap_int, exc_wrap_pointer
@ -62,12 +63,12 @@ def _shapes(image, mask, connectivity, transform):
cdef ShapeIterator shape_iter = None
cdef int fieldtp
is_float = np.dtype(image.dtype).kind == "f"
is_float = _getnpdtype(image.dtype).kind == "f"
fieldtp = 2 if is_float else 0
valid_dtypes = ('int16', 'int32', 'uint8', 'uint16', 'float32')
if np.dtype(image.dtype).name not in valid_dtypes:
if _getnpdtype(image.dtype).name not in valid_dtypes:
raise ValueError("image dtype must be one of: {0}".format(
', '.join(valid_dtypes)))
@ -89,7 +90,7 @@ def _shapes(image, mask, connectivity, transform):
if mask.shape != image.shape:
raise ValueError("Mask must have same shape as image")
if np.dtype(mask.dtype).name not in ('bool', 'uint8'):
if _getnpdtype(mask.dtype).name not in ('bool', 'uint8'):
raise ValueError("Mask must be dtype rasterio.bool_ or "
"rasterio.uint8")
@ -184,7 +185,7 @@ def _sieve(image, size, out, mask, connectivity):
valid_dtypes = ('int16', 'int32', 'uint8', 'uint16')
if np.dtype(image.dtype).name not in valid_dtypes:
if _getnpdtype(image.dtype).name not in valid_dtypes:
valid_types_str = ', '.join(('rasterio.{0}'.format(t) for t
in valid_dtypes))
raise ValueError(
@ -203,7 +204,7 @@ def _sieve(image, size, out, mask, connectivity):
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:
if _getnpdtype(image.dtype).name != _getnpdtype(out.dtype).name:
raise ValueError('out raster must match dtype of image')
try:
@ -231,7 +232,7 @@ def _sieve(image, size, out, mask, connectivity):
if mask.shape != image.shape:
raise ValueError("Mask must have same shape as image")
if np.dtype(mask.dtype) not in ('bool', 'uint8'):
if _getnpdtype(mask.dtype) not in ('bool', 'uint8'):
raise ValueError("Mask must be dtype rasterio.bool_ or "
"rasterio.uint8")

View File

@ -26,7 +26,7 @@ from rasterio.errors import (
NotGeoreferencedWarning, NodataShadowWarning, WindowError,
UnsupportedOperation, OverviewCreationError, RasterBlockError, InvalidArrayError
)
from rasterio.dtypes import is_ndarray
from rasterio.dtypes import is_ndarray, _is_complex_int, _getnpdtype
from rasterio.sample import sample_gen
from rasterio.transform import Affine
from rasterio.path import parse_path, UnparsedPath
@ -300,7 +300,7 @@ cdef bint in_dtype_range(value, dtype):
105: np.iinfo,
117: np.iinfo}
key = np.dtype(dtype).kind
key = _getnpdtype(dtype).kind
if np.isnan(value):
return key in ('c', 'f', 99, 102)
@ -437,12 +437,14 @@ cdef class DatasetReaderBase(DatasetBase):
check_dtypes = set()
nodatavals = []
# Check each index before processing 3D array
for bidx in indexes:
if bidx not in self.indexes:
raise IndexError("band index {} out of range (not in {})".format(bidx, self.indexes))
idx = self.indexes.index(bidx)
idx = self.indexes.index(bidx)
dtype = self.dtypes[idx]
check_dtypes.add(dtype)
@ -450,16 +452,18 @@ cdef class DatasetReaderBase(DatasetBase):
log.debug("Output nodata value read from file: %r", ndv)
if ndv is not None:
kind = np.dtype(dtype).kind
if chr(kind) in "iu":
if ndv is not None and not _is_complex_int(dtype):
kind = _getnpdtype(dtype).kind
if kind in "iu":
info = np.iinfo(dtype)
dt_min, dt_max = info.min, info.max
elif chr(kind) in "cf":
elif kind in "cf":
info = np.finfo(dtype)
dt_min, dt_max = info.min, info.max
else:
dt_min, dt_max = False, True
if ndv < dt_min:
ndv = dt_min
elif ndv > dt_max:
@ -480,6 +484,9 @@ cdef class DatasetReaderBase(DatasetBase):
if out_dtype is not None:
dtype = out_dtype
# Ensure we have a numpy dtype.
dtype = _getnpdtype(dtype)
# Get the natural shape of the read window, boundless or not.
# The window can have float values. In this case, we round up
# when computing the shape.
@ -553,8 +560,7 @@ cdef class DatasetReaderBase(DatasetBase):
log.debug("Jump straight to _read()")
log.debug("Window: %r", window)
out = self._read(indexes, out, window, dtype,
resampling=resampling)
out = self._read(indexes, out, window, dtype, resampling=resampling)
if masked or fill_value is not None:
if all_valid:
@ -751,7 +757,7 @@ cdef class DatasetReaderBase(DatasetBase):
out = np.zeros(out_shape, 'uint8')
if out is not None:
if out.dtype != np.dtype(dtype):
if out.dtype != _getnpdtype(dtype):
raise ValueError(
"the out array's dtype '%s' does not match '%s'"
% (out.dtype, dtype))
@ -820,8 +826,7 @@ cdef class DatasetReaderBase(DatasetBase):
return out
def _read(self, indexes, out, window, dtype, masks=False,
resampling=Resampling.nearest):
def _read(self, indexes, out, window, dtype, masks=False, resampling=Resampling.nearest):
"""Read raster bands as a multidimensional array
If `indexes` is a list, the result is a 3D array, but
@ -1252,19 +1257,21 @@ cdef class DatasetWriterBase(DatasetReaderBase):
try:
width = int(width)
height = int(height)
except:
except TypeError:
raise TypeError("Integer width and height are required.")
try:
count = int(count)
except:
except TypeError:
raise TypeError("Integer band count is required.")
try:
assert dtype is not None
_ = np.dtype(dtype)
except:
raise TypeError("A valid dtype is required.")
self._init_dtype = np.dtype(dtype).name
if _is_complex_int(dtype):
self._init_dtype = dtype
else:
try:
assert dtype is not None
self._init_dtype = _getnpdtype(dtype).name
except Exception:
raise TypeError("A valid dtype is required.")
# Make and store a GDAL dataset handle.
filename = path.name
@ -1343,7 +1350,9 @@ cdef class DatasetWriterBase(DatasetReaderBase):
if nodata is not None:
if not in_dtype_range(nodata, dtype):
if _is_complex_int(dtype):
pass
elif not in_dtype_range(nodata, dtype):
raise ValueError(
"Given nodata value, %s, is beyond the valid "
"range of its data type, %s." % (
@ -1561,10 +1570,7 @@ cdef class DatasetWriterBase(DatasetReaderBase):
else: # unique dtype; normal case
dtype = check_dtypes.pop()
if arr is not None and arr.dtype != dtype:
raise ValueError(
"the array's dtype '%s' does not match "
"the file's dtype '%s'" % (arr.dtype, dtype))
dtype = _getnpdtype(dtype)
# Require C-continguous arrays (see #108).
arr = np.require(arr, dtype=dtype, requirements='C')
@ -1769,7 +1775,7 @@ cdef class DatasetWriterBase(DatasetReaderBase):
GDALFillRaster(mask, 255, 0)
elif mask_array is False:
GDALFillRaster(mask, 0, 0)
elif mask_array.dtype == np.bool:
elif mask_array.dtype == bool:
array = 255 * mask_array.astype(np.uint8)
io_band(mask, 1, xoff, yoff, width, height, array)
else:
@ -2163,13 +2169,14 @@ cdef class BufferedDatasetWriterBase(DatasetWriterBase):
count = int(count)
except:
raise TypeError("Integer band count is required.")
try:
assert dtype is not None
_ = np.dtype(dtype)
except:
_ = _getnpdtype(dtype)
except Exception:
raise TypeError("A valid dtype is required.")
self._init_dtype = np.dtype(dtype).name
self._init_dtype = _getnpdtype(dtype).name
self.name = path.name
self.mode = mode

View File

@ -16,7 +16,7 @@ from rasterio.enums import Resampling
cimport numpy as np
from rasterio._err cimport exc_wrap_int, exc_wrap_pointer
from rasterio.errors import GDALOptionNotImplementedError
from rasterio.errors import GDALOptionNotImplementedError, DatasetIOShapeError
cdef GDALDatasetH open_dataset(
@ -114,6 +114,9 @@ cdef int io_multi_band(
cdef int xsize = <int>width
cdef int ysize = <int>height
if len(indexes) != data.shape[0]:
raise DatasetIOShapeError("Dataset indexes and destination buffer are mismatched")
bandmap = <int *>CPLMalloc(count*sizeof(int))
for i in range(count):
bandmap[i] = <int>indexes[i]
@ -162,6 +165,9 @@ cdef int io_multi_mask(
cdef int xsize = <int>width
cdef int ysize = <int>height
if len(indexes) != data.shape[0]:
raise DatasetIOShapeError("Dataset indexes and destination buffer are mismatched")
for i in range(count):
j = <int>indexes[i]
band = GDALGetRasterBand(hds, j)

View File

@ -43,6 +43,17 @@ from rasterio._features cimport GeomBuilder, OGRGeomBuilder
log = logging.getLogger(__name__)
# Gauss (7) is not supported for warp
SUPPORTED_RESAMPLING = [r for r in Resampling if r.value < 7]
GDAL2_RESAMPLING = [r for r in Resampling if r.value > 7 and r.value <= 12]
if GDALVersion.runtime().at_least('2.0'):
SUPPORTED_RESAMPLING.extend(GDAL2_RESAMPLING)
# sum supported since GDAL 3.1
if GDALVersion.runtime().at_least('3.1'):
SUPPORTED_RESAMPLING.append(Resampling.sum)
# rms supported since GDAL 3.3
if GDALVersion.runtime().at_least('3.3'):
SUPPORTED_RESAMPLING.append(Resampling.rms)
def recursive_round(val, precision):
"""Recursively round coordinates."""
@ -774,6 +785,16 @@ cdef class WarpedVRTReaderBase(DatasetReaderBase):
if src_dataset.mode != "r":
warnings.warn("Source dataset should be opened in read-only mode. Use of datasets opened in modes other than 'r' will be disallowed in a future version.", RasterioDeprecationWarning, stacklevel=2)
# Guard against invalid or unsupported resampling algorithms.
try:
if resampling == 7:
raise ValueError("Gauss resampling is not supported")
Resampling(resampling)
except ValueError:
raise ValueError(
"resampling must be one of: {0}".format(", ".join(['Resampling.{0}'.format(r.name) for r in SUPPORTED_RESAMPLING])))
self.mode = 'r'
self.options = {}
self._count = 0

View File

@ -1,11 +1,9 @@
"""Bounding box tuple, and disjoint operator."""
from collections import namedtuple, OrderedDict
from collections import namedtuple
_BoundingBox = namedtuple('BoundingBox', ('left', 'bottom', 'right', 'top'))
class BoundingBox(_BoundingBox):
BoundingBox = namedtuple('BoundingBox', ('left', 'bottom', 'right', 'top'))
BoundingBox.__doc__ = \
"""Bounding box named tuple, defining extent in cartesian coordinates.
.. code::
@ -24,10 +22,6 @@ class BoundingBox(_BoundingBox):
Top coordinate
"""
def _asdict(self):
return OrderedDict(zip(self._fields, self))
def disjoint_bounds(bounds1, bounds2):
"""Compare two bounds and determine if they are disjoint.
@ -43,6 +37,7 @@ def disjoint_bounds(bounds1, bounds2):
boolean
``True`` if bounds are disjoint,
``False`` if bounds overlap
"""
bounds1_north_up = bounds1[3] > bounds1[1]
bounds2_north_up = bounds2[3] > bounds2[1]

View File

@ -4,10 +4,6 @@ Since 0.13 we are not importing numpy here and data types are strings.
Happily strings can be used throughout Numpy and so existing code will
not break.
Within Rasterio, to test data types, we use Numpy's dtype() factory to
do something like this:
if np.dtype(destination.dtype) == np.dtype(rasterio.uint8): ...
"""
bool_ = 'bool'
@ -23,26 +19,29 @@ complex_ = 'complex'
complex64 = 'complex64'
complex128 = 'complex128'
# Not supported:
# GDT_CInt16 = 8, GDT_CInt32 = 9, GDT_CFloat32 = 10, GDT_CFloat64 = 11
complex_int16 = "complex_int16"
dtype_fwd = {
0: None, # GDT_Unknown
1: ubyte, # GDT_Byte
2: uint16, # GDT_UInt16
3: int16, # GDT_Int16
4: uint32, # GDT_UInt32
5: int32, # GDT_Int32
6: float32, # GDT_Float32
7: float64, # GDT_Float64
8: complex_, # GDT_CInt16
9: complex_, # GDT_CInt32
10: complex64, # GDT_CFloat32
11: complex128} # GDT_CFloat64
0: None, # GDT_Unknown
1: ubyte, # GDT_Byte
2: uint16, # GDT_UInt16
3: int16, # GDT_Int16
4: uint32, # GDT_UInt32
5: int32, # GDT_Int32
6: float32, # GDT_Float32
7: float64, # GDT_Float64
8: complex_int16, # GDT_CInt16
9: complex64, # GDT_CInt32
10: complex64, # GDT_CFloat32
11: complex128, # GDT_CFloat64
}
dtype_rev = dict((v, k) for k, v in dtype_fwd.items())
dtype_rev['uint8'] = 1
dtype_rev['int8'] = 1
dtype_rev["uint8"] = 1
dtype_rev["int8"] = 1
dtype_rev["complex"] = 11
dtype_rev["complex_int16"] = 8
typename_fwd = {
0: 'Unknown',
@ -154,7 +153,7 @@ def can_cast_dtype(values, dtype):
if not is_ndarray(values):
values = np.array(values)
if values.dtype.name == np.dtype(dtype).name:
if values.dtype.name == _getnpdtype(dtype).name:
return True
elif values.dtype.kind == 'f':
@ -185,3 +184,15 @@ def validate_dtype(values, valid_dtypes):
return (values.dtype.name in valid_dtypes or
get_minimum_dtype(values) in valid_dtypes)
def _is_complex_int(dtype):
return isinstance(dtype, str) and dtype.startswith("complex_int")
def _getnpdtype(dtype):
import numpy as np
if _is_complex_int(dtype):
return np.dtype("complex64")
else:
return np.dtype(dtype)

View File

@ -87,9 +87,9 @@ class GDALOptionNotImplementedError(RasterioError):
by GDAL 1.x.
"""
class GDALVersionError(RasterioError):
"""Raised if the runtime version of GDAL does not meet the required
version of GDAL."""
"""Raised if the runtime version of GDAL does not meet the required version of GDAL."""
class WindowEvaluationError(ValueError):
@ -142,3 +142,7 @@ class TransformError(RasterioError):
class WarpedVRTError(RasterioError):
"""Raised when WarpedVRT can't be initialized"""
class DatasetIOShapeError(RasterioError):
"""Raised when data buffer shape is a mismatch when reading and writing"""

View File

@ -1,6 +1,5 @@
"""Functions for working with features in a raster dataset."""
import logging
import math
import os
@ -8,11 +7,11 @@ import warnings
import numpy as np
import rasterio._loading
with rasterio._loading.add_gdal_dll_directories():
import rasterio
from rasterio.dtypes import validate_dtype, can_cast_dtype, get_minimum_dtype
from rasterio.dtypes import validate_dtype, can_cast_dtype, get_minimum_dtype, getnpdtype
from rasterio.enums import MergeAlg
from rasterio.env import ensure_env
from rasterio.errors import ShapeSkipWarning
@ -279,7 +278,7 @@ def rasterize(
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:
if dtype is not None and getnpdtype(dtype).name not in valid_dtypes:
raise ValueError(format_invalid_dtype('dtype'))
valid_shapes = []
@ -335,7 +334,7 @@ def rasterize(
raise ValueError(format_cast_error('shape values', dtype))
if out is not None:
if np.dtype(out.dtype).name not in valid_dtypes:
if getnpdtype(out.dtype).name not in valid_dtypes:
raise ValueError(format_invalid_dtype('out'))
if not can_cast_dtype(shape_values, out.dtype):

View File

@ -9,8 +9,10 @@ import warnings
import numpy as np
import rasterio._loading
with rasterio._loading.add_gdal_dll_directories():
import rasterio
from rasterio.coords import disjoint_bounds
from rasterio.enums import Resampling
from rasterio import windows
from rasterio.transform import Affine
@ -23,13 +25,13 @@ def copy_first(merged_data, new_data, merged_mask, new_mask, **kwargs):
mask = np.empty_like(merged_mask, dtype="bool")
np.logical_not(new_mask, out=mask)
np.logical_and(merged_mask, mask, out=mask)
np.copyto(merged_data, new_data, where=mask)
np.copyto(merged_data, new_data, where=mask, casting="unsafe")
def copy_last(merged_data, new_data, merged_mask, new_mask, **kwargs):
mask = np.empty_like(merged_mask, dtype="bool")
np.logical_not(new_mask, out=mask)
np.copyto(merged_data, new_data, where=mask)
np.copyto(merged_data, new_data, where=mask, casting="unsafe")
def copy_min(merged_data, new_data, merged_mask, new_mask, **kwargs):
@ -39,7 +41,7 @@ def copy_min(merged_data, new_data, merged_mask, new_mask, **kwargs):
np.minimum(merged_data, new_data, out=merged_data, where=mask)
np.logical_not(new_mask, out=mask)
np.logical_and(merged_mask, mask, out=mask)
np.copyto(merged_data, new_data, where=mask)
np.copyto(merged_data, new_data, where=mask, casting="unsafe")
def copy_max(merged_data, new_data, merged_mask, new_mask, **kwargs):
@ -49,7 +51,7 @@ def copy_max(merged_data, new_data, merged_mask, new_mask, **kwargs):
np.maximum(merged_data, new_data, out=merged_data, where=mask)
np.logical_not(new_mask, out=mask)
np.logical_and(merged_mask, mask, out=mask)
np.copyto(merged_data, new_data, where=mask)
np.copyto(merged_data, new_data, where=mask, casting="unsafe")
MERGE_METHODS = {
@ -292,6 +294,10 @@ def merge(
# This approach uses the maximum amount of memory to solve the
# problem. Making it more efficient is a TODO.
if disjoint_bounds((dst_w, dst_s, dst_e, dst_n), src.bounds):
logger.debug("Skipping source: src=%r, window=%r", src)
continue
# 1. Compute spatial intersection of destination and source
src_w, src_s, src_e, src_n = src.bounds
@ -304,7 +310,6 @@ def merge(
src_window = windows.from_bounds(
int_w, int_s, int_e, int_n, src.transform, precision=precision
)
logger.debug("Src %s window: %r", src.name, src_window)
# 3. Compute the destination window
dst_window = windows.from_bounds(

View File

@ -190,13 +190,13 @@ def edit(ctx, input, bidx, nodata, unset_nodata, crs, unset_crs, transform,
tags = allmd['tags']
colorinterp = allmd['colorinterp']
if unset_nodata and nodata is not options.IgnoreOption:
if unset_nodata and nodata is not None:
raise click.BadParameter(
"--unset-nodata and --nodata cannot be used together.")
"--unset-nodata and --nodata cannot be used together."
)
if unset_crs and crs:
raise click.BadParameter(
"--unset-crs and --crs cannot be used together.")
raise click.BadParameter("--unset-crs and --crs cannot be used together.")
if unset_nodata:
# Setting nodata to None will raise NotImplementedError
@ -207,17 +207,18 @@ def edit(ctx, input, bidx, nodata, unset_nodata, crs, unset_crs, transform,
except NotImplementedError as exc: # pragma: no cover
raise click.ClickException(str(exc))
elif nodata is not options.IgnoreOption:
elif nodata is not None:
dtype = dst.dtypes[0]
if nodata is not None and not in_dtype_range(nodata, dtype):
raise click.BadParameter(
"outside the range of the file's "
"data type (%s)." % dtype,
param=nodata, param_hint='nodata')
"outside the range of the file's data type (%s)." % dtype,
param=nodata,
param_hint="nodata",
)
dst.nodata = nodata
if unset_crs:
dst.crs = None # CRS()
dst.crs = None
elif crs:
dst.crs = crs

View File

@ -51,26 +51,12 @@ import click
import rasterio
import rasterio.shutil
from rasterio.path import parse_path, ParsedPath, UnparsedPath
from rasterio.path import parse_path, UnparsedPath
logger = logging.getLogger(__name__)
class IgnoreOptionMarker:
"""A marker for an option that is to be ignored.
For use in the case where `None` is a meaningful option value,
such as for nodata where `None` means "unset the nodata value."
"""
def __repr__(self):
return 'IgnoreOption()'
IgnoreOption = IgnoreOptionMarker()
def _cb_key_val(ctx, param, value):
"""
@ -180,29 +166,36 @@ def like_handler(ctx, param, value):
def nodata_handler(ctx, param, value):
"""Return a float or None"""
if value is None or value is IgnoreOption:
return value
elif value.lower() in ['null', 'nil', 'none', 'nada']:
if value is None or value.lower() in ["null", "nil", "none", "nada"]:
return None
else:
try:
return float(value)
except (TypeError, ValueError):
raise click.BadParameter(
"{!r} is not a number".format(value),
param=param, param_hint='nodata')
"{!r} is not a number".format(value), param=param, param_hint="nodata"
)
def edit_nodata_handler(ctx, param, value):
"""Get nodata value from a template file or command line.
Expected values are 'like', 'null', a numeric value, 'nan', or
IgnoreOption. Anything else should raise BadParameter.
Expected values are 'like', 'null', a numeric value, 'nan', or None.
Returns
-------
float or None
Raises
------
click.BadParameter
"""
if value == 'like' or value is IgnoreOption:
if value == "like" or value is None:
retval = from_like_context(ctx, param, value)
if retval is not None:
return retval
return nodata_handler(ctx, param, value)
@ -335,12 +328,20 @@ overwrite_opt = click.option(
help="Always overwrite an existing output file.")
nodata_opt = click.option(
'--nodata', callback=nodata_handler, default=None,
metavar='NUMBER|nan', help="Set a Nodata value.")
"--nodata",
callback=nodata_handler,
default=None,
metavar="NUMBER|nan",
help="Set a Nodata value.",
)
edit_nodata_opt = click.option(
'--nodata', callback=edit_nodata_handler, default=IgnoreOption,
metavar='NUMBER|nan|null', help="Modify the Nodata value.")
"--nodata",
callback=edit_nodata_handler,
default=None,
metavar="NUMBER|nan|null",
help="Modify the Nodata value.",
)
like_opt = click.option(
'--like',

View File

@ -189,9 +189,14 @@ def union(*windows):
Window
"""
stacked = np.dstack([toranges(w) for w in windows])
return Window.from_slices(
(stacked[0, 0].min(), stacked[0, 1].max()),
(stacked[1, 0].min(), stacked[1, 1]. max()))
row_start, row_stop = stacked[0, 0].min(), stacked[0, 1].max()
col_start, col_stop = stacked[1, 0].min(), stacked[1, 1].max()
return Window(
col_off=col_start,
row_off=row_start,
width=col_stop - col_start,
height=row_stop - row_start,
)
@iter_args
@ -213,9 +218,14 @@ def intersection(*windows):
raise WindowError("windows do not intersect")
stacked = np.dstack([toranges(w) for w in windows])
return Window.from_slices(
(stacked[0, 0].max(), stacked[0, 1].min()),
(stacked[1, 0].max(), stacked[1, 1]. min()))
row_start, row_stop = stacked[0, 0].max(), stacked[0, 1].min()
col_start, col_stop = stacked[1, 0].max(), stacked[1, 1].min()
return Window(
col_off=col_start,
row_off=row_start,
width=col_stop - col_start,
height=row_stop - row_start,
)
@iter_args

View File

@ -255,7 +255,7 @@ inst_reqs = [
"affine",
"attrs",
"certifi",
"click>=4.0,<8",
"click>=4.0",
"cligj>=0.5",
"numpy",
"snuggs>=1.4.1",

View File

@ -1,4 +1,5 @@
import logging
import subprocess
import sys
import uuid
@ -8,9 +9,6 @@ import pytest
import rasterio
logging.basicConfig(stream=sys.stderr, level=logging.DEBUG)
@pytest.fixture(scope='function')
def tempfile():
"""A temporary filename in the GDAL '/vsimem' filesystem"""
@ -23,10 +21,8 @@ def complex_image(height, width, dtype):
[complex(x, x) for x in range(height * width)],
dtype=dtype).reshape(height, width)
dtypes = ['complex', 'complex64', 'complex128']
@pytest.mark.parametrize("dtype", dtypes)
@pytest.mark.parametrize("dtype", ["complex", "complex64", "complex128"])
@pytest.mark.parametrize("height,width", [(20, 20)])
def test_read_array(tempfile, dtype, height, width):
"""_io functions read and write arrays correctly"""
@ -58,3 +54,45 @@ def test_complex_nodata(tmpdir):
with rasterio.open(tempfile) as dst:
assert dst.nodata == 0
@pytest.mark.gdalbin
def test_complex_int16(tmpdir):
"""A cint16 dataset can be created"""
import numpy as np
import rasterio
from rasterio.transform import Affine
x = np.linspace(-4.0, 4.0, 240)
y = np.linspace(-3.0, 3.0, 180)
X, Y = np.meshgrid(x, y)
Z1 = np.ones_like(X) + 1j
res = (x[-1] - x[0]) / 240.0
transform1 = Affine.translation(x[0] - res / 2, y[-1] - res / 2) * Affine.scale(
res, -res
)
tempfile = str(tmpdir.join("test.tif"))
with rasterio.open(
tempfile,
"w",
driver="GTiff",
height=Z1.shape[0],
width=Z1.shape[1],
nodata=0,
count=1,
dtype="complex_int16",
crs="+proj=latlong",
transform=transform1,
) as dst:
dst.write(Z1, 1)
assert "Type=CInt16" in subprocess.check_output(["gdalinfo", tempfile]).decode(
"utf-8"
)
with rasterio.open(tempfile) as dst:
data = dst.read()
assert data.dtype == np.complex64

View File

@ -3,10 +3,26 @@ import pytest
import rasterio
from rasterio import (
ubyte, uint8, uint16, uint32, int16, int32, float32, float64, complex_)
ubyte,
uint8,
uint16,
uint32,
int16,
int32,
float32,
float64,
complex_,
complex_int16,
)
from rasterio.dtypes import (
_gdal_typename, is_ndarray, check_dtype, get_minimum_dtype, can_cast_dtype,
validate_dtype
_gdal_typename,
is_ndarray,
check_dtype,
get_minimum_dtype,
can_cast_dtype,
validate_dtype,
_is_complex_int,
_getnpdtype,
)
@ -28,10 +44,19 @@ def test_check_dtype_invalid():
assert not check_dtype('foo')
def test_gdal_name():
assert _gdal_typename(ubyte) == 'Byte'
assert _gdal_typename(np.uint8) == 'Byte'
assert _gdal_typename(np.uint16) == 'UInt16'
@pytest.mark.parametrize(
("dtype", "name"),
[
(ubyte, "Byte"),
(np.uint8, "Byte"),
(np.uint16, "UInt16"),
("uint8", "Byte"),
("complex_int16", "CInt16"),
(complex_int16, "CInt16"),
],
)
def test_gdal_name(dtype, name):
assert _gdal_typename(dtype) == name
def test_get_minimum_dtype():
@ -46,8 +71,8 @@ def test_get_minimum_dtype():
assert get_minimum_dtype(np.array([0, 1], dtype=np.uint)) == uint8
assert get_minimum_dtype(np.array([0, 1000], dtype=np.uint)) == uint16
assert get_minimum_dtype(np.array([0, 100000], dtype=np.uint)) == uint32
assert get_minimum_dtype(np.array([-1, 0, 1], dtype=np.int)) == int16
assert get_minimum_dtype(np.array([-1, 0, 100000], dtype=np.int)) == int32
assert get_minimum_dtype(np.array([-1, 0, 1], dtype=int)) == int16
assert get_minimum_dtype(np.array([-1, 0, 100000], dtype=int)) == int32
assert get_minimum_dtype(np.array([-1.5, 0, 1.5], dtype=np.float64)) == float32
@ -89,3 +114,17 @@ def test_complex(tmpdir):
arr2 = src.read(1)
assert np.array_equal(arr1, arr2)
def test_is_complex_int():
assert _is_complex_int("complex_int16")
def test_not_is_complex_int():
assert not _is_complex_int("complex")
def test_get_npdtype():
npdtype = _getnpdtype("complex_int16")
assert npdtype == np.complex64
assert npdtype.kind == "c"

View File

@ -59,3 +59,10 @@ def test_issue2163():
data = src.read()
result, transform = merge([src])
assert numpy.allclose(data, result)
def test_unsafe_casting():
"""Demonstrate fix for issue 2179"""
with rasterio.open("tests/data/float_raster_with_nodata.tif") as src:
result, transform = merge([src], dtype="uint8", nodata=0.0)
assert not result.any() # this is why it's called "unsafe".

View File

@ -1,23 +1,18 @@
from hashlib import md5
import logging
import sys
import unittest
import numpy as np
import pytest
import rasterio
logging.basicConfig(stream=sys.stderr, level=logging.DEBUG)
from rasterio.errors import DatasetIOShapeError
# Find out if we've got HDF support (needed below).
try:
with rasterio.open('tests/data/no_band.h5') as s:
pass
has_hdf = True
except:
except Exception:
has_hdf = False
@ -93,13 +88,7 @@ class ReaderContextTest(unittest.TestCase):
def test_read_out_dtype_fail(self):
with rasterio.open('tests/data/RGB.byte.tif') as s:
a = np.zeros((718, 791), dtype=rasterio.float32)
try:
s.read(1, a)
except ValueError as e:
assert ("the array's dtype 'float32' does not match the "
"file's dtype") in str(e)
except:
assert "failed to catch exception" is False
s.read(1, a)
def test_read_basic(self):
with rasterio.open('tests/data/shade.tif') as s:
@ -315,7 +304,7 @@ def test_out_shape_exceptions(path_rgb_byte_tif):
out_shape = (src.count, src.height, src.width)
reader(out=out, out_shape=out_shape)
with pytest.raises(ValueError):
with pytest.raises(Exception):
out_shape = (5, src.height, src.width)
reader(1, out_shape=out_shape)
@ -325,3 +314,10 @@ def test_out_shape_implicit(path_rgb_byte_tif):
with rasterio.open(path_rgb_byte_tif) as src:
out = src.read(indexes=(1, 2), out_shape=src.shape)
assert out.shape == (2,) + src.shape
def test_out_shape_no_segfault(path_rgb_byte_tif):
"""Prevent segfault as reported in 2189"""
with rasterio.open(path_rgb_byte_tif) as src:
with pytest.raises(DatasetIOShapeError):
src.read(out_shape=(2, src.height, src.width))

View File

@ -65,7 +65,7 @@ def test_delete_nodata(data, runner):
"""Delete a dataset's nodata value"""
inputfile = str(data.join('RGB.byte.tif'))
result = runner.invoke(
main_group, ['edit-info', inputfile, '--unset-nodata'])
main_group, ['edit-info', inputfile, '--unset-nodata'], catch_exceptions=False)
assert result.exit_code == 0

View File

@ -9,8 +9,13 @@ import pytest
from rasterio.enums import ColorInterp
from rasterio.rio.options import (
IgnoreOption, bounds_handler, file_in_handler, like_handler,
edit_nodata_handler, nodata_handler, _cb_key_val)
bounds_handler,
file_in_handler,
like_handler,
edit_nodata_handler,
nodata_handler,
_cb_key_val,
)
class MockContext:
@ -186,14 +191,14 @@ def test_edit_nodata_callback_like(data):
def test_edit_nodata_callback_all_like(data):
ctx = MockContext()
ctx.obj['like'] = {'nodata': 0.0}
ctx.obj['all_like'] = True
assert edit_nodata_handler(ctx, MockOption('nodata'), IgnoreOption) == 0.0
ctx.obj["like"] = {"nodata": 0.0}
ctx.obj["all_like"] = True
assert edit_nodata_handler(ctx, MockOption("nodata"), None) == 0.0
def test_edit_nodata_callback_ignore(data):
ctx = MockContext()
assert edit_nodata_handler(ctx, MockOption('nodata'), IgnoreOption) is IgnoreOption
assert edit_nodata_handler(ctx, MockOption("nodata"), None) is None
def test_edit_nodata_callback_none(data):

View File

@ -1,7 +1,6 @@
"""rasterio.warp module tests"""
import json
import sys
import logging
from affine import Affine
@ -15,7 +14,6 @@ from rasterio.crs import CRS
from rasterio.enums import Resampling
from rasterio.env import GDALVersion
from rasterio.errors import (
GDALBehaviorChangeException,
CRSError,
GDALVersionError,
TransformError,
@ -1136,6 +1134,7 @@ def test_reproject_resampling(path_rgb_byte_tif, method):
assert np.count_nonzero(out) in expected[method]
@pytest.mark.parametrize("test3d,count_nonzero", [(True, 1309625), (False, 437686)])
def test_reproject_array_interface(test3d, count_nonzero, path_rgb_byte_tif):
class DataArray:
@ -1333,11 +1332,12 @@ def test_resample_default_invert_proj(method):
assert out.mean() > 0
@pytest.mark.xfail(reason="Projection extents have changed with PROJ 8")
def test_target_aligned_pixels():
"""Issue 853 has been resolved"""
with rasterio.open("tests/data/world.rgb.tif") as src:
source = src.read(1)
profile = src.profile.copy()
profile = src.profile
dst_crs = "EPSG:3857"
@ -1348,7 +1348,7 @@ def test_target_aligned_pixels():
)
dst_affine, dst_width, dst_height = aligned_target(
dst_affine, dst_width, dst_height, 100000.0
dst_affine, dst_width, dst_height, 10000.0
)
profile["height"] = dst_height
@ -1366,7 +1366,7 @@ def test_target_aligned_pixels():
resampling=Resampling.nearest,
)
# Check that there is no black borders
# Check that there are no black borders
assert out[:, 0].all()
assert out[:, -1].all()
assert out[0, :].all()
@ -1694,26 +1694,6 @@ def test_issue_1446_b():
assert all([-350 < x < -150 for x, y in transformed_geoms[183519]["coordinates"]])
def test_issue_1076():
"""Confirm fix of #1076"""
arr = (np.random.random((20, 30)) * 100).astype('int32')
fill_value = 42
newarr = np.full((200, 300), fill_value=fill_value, dtype='int32')
src_crs = CRS.from_epsg(32632)
src_transform = Affine(600.0, 0.0, 399960.0, 0.0, -600.0, 6100020.0)
dst_transform = Affine(60.0, 0.0, 399960.0, 0.0, -60.0, 6100020.0)
reproject(arr, newarr,
src_transform=src_transform,
dst_transform=dst_transform,
src_crs=src_crs,
dst_crs=src_crs,
resample=Resampling.nearest)
assert not (newarr == fill_value).all()
def test_reproject_init_dest_nodata():
"""No pixels should transfer over"""
crs = CRS.from_epsg(4326)
@ -1785,16 +1765,16 @@ def test_reproject_rpcs_with_transformer_options(caplog):
transform = Affine(0.001953396267361111, 0.0, -124.00013888888888, 0.0, -0.001953396267361111, 50.000138888888884)
with mem.open(
driver="GTiff",
width=1024,
height=1024,
width=1024,
height=1024,
count=1,
transform=transform,
dtype='int16',
crs=crs
dtype="int16",
crs=crs,
) as dem:
# we flush dem dataset before letting GDAL read from vsimem
dem.write_band(1, 500 * np.ones((1024, 1024), dtype='int16'))
out = np.zeros(
(3, src.profile["width"], src.profile["height"]), dtype=np.uint8
)
@ -1810,7 +1790,7 @@ def test_reproject_rpcs_with_transformer_options(caplog):
resampling=Resampling.nearest,
RPC_DEM=dem.name,
)
)
caplog.set_level(logging.INFO)
reproject(
rasterio.band(src, src.indexes),
@ -1820,8 +1800,8 @@ def test_reproject_rpcs_with_transformer_options(caplog):
dst_crs="EPSG:3857",
resampling=Resampling.nearest,
)
)
assert not out.all()
assert not out2.all()
assert "RPC_DEM" in caplog.text
@ -1853,10 +1833,10 @@ def test_warp_gcps_compute_dst_transform_automatically_array():
assert not out[:, -1, -1].any()
assert not out[:, -1, 0].any()
def test_warp_gcps_compute_dst_transform_automatically_reader(tmpdir):
"""Ensure we don't raise an exception when gcps passed without dst_transform, for a source dataset"""
tiffname = str(tmpdir.join('test.tif'))
arr = np.ones((3, 800, 800), dtype=np.uint8) * 255
src_gcps = [
GroundControlPoint(row=0, col=0, x=156113, y=2818720, z=0),
GroundControlPoint(row=0, col=800, x=338353, y=2785790, z=0),
@ -1867,7 +1847,7 @@ def test_warp_gcps_compute_dst_transform_automatically_reader(tmpdir):
with rasterio.open(tiffname, mode='w', height=800, width=800, count=3, dtype=np.uint8) as source:
source.gcps = (src_gcps, CRS.from_epsg(32618))
with rasterio.open(tiffname) as source:
reproject(
rasterio.band(source, source.indexes),
@ -1875,13 +1855,14 @@ def test_warp_gcps_compute_dst_transform_automatically_reader(tmpdir):
dst_crs="EPSG:32618",
resampling=Resampling.nearest
)
assert not out.all()
assert not out[:, 0, 0].any()
assert not out[:, 0, -1].any()
assert not out[:, -1, -1].any()
assert not out[:, -1, 0].any()
def test_reproject_rpcs_exact_transformer(caplog):
"""Reproject using rational polynomial coefficients and DEM, requiring that
we don't try to make an approximate transformer.
@ -1892,16 +1873,16 @@ def test_reproject_rpcs_exact_transformer(caplog):
transform = Affine(0.001953396267361111, 0.0, -124.00013888888888, 0.0, -0.001953396267361111, 50.000138888888884)
with mem.open(
driver="GTiff",
width=1024,
height=1024,
width=1024,
height=1024,
count=1,
transform=transform,
dtype='int16',
crs=crs
dtype="int16",
crs=crs,
) as dem:
# we flush dem dataset before letting GDAL read from vsimem
dem.write_band(1, 500 * np.ones((1024, 1024), dtype='int16'))
out = np.zeros(
(3, src.profile["width"], src.profile["height"]), dtype=np.uint8
)
@ -1915,10 +1896,10 @@ def test_reproject_rpcs_exact_transformer(caplog):
dst_crs="EPSG:3857",
resampling=Resampling.nearest,
RPC_DEM=dem.name,
)
)
assert "Created exact transformer" in caplog.text
def test_reproject_rpcs_approx_transformer(caplog):
"""Reproject using rational polynomial coefficients without a DEM, for which it's
@ -1937,6 +1918,6 @@ def test_reproject_rpcs_approx_transformer(caplog):
rpcs=src_rpcs,
dst_crs="EPSG:3857",
resampling=Resampling.nearest,
)
)
assert "Created approximate transformer" in caplog.text
assert "Created approximate transformer" in caplog.text

View File

@ -608,3 +608,10 @@ def test_issue2086():
with rasterio.open("tests/data/white-gemini-iv.vrt") as src:
with WarpedVRT(src, crs=DST_CRS) as vrt:
assert vrt.shape == (1031, 1146)
def test_gauss_no(path_rgb_byte_tif):
"""Guard against the issue reported in #2190"""
with rasterio.open(path_rgb_byte_tif) as src:
with pytest.raises(Exception):
WarpedVRT(src, resampling=Resampling.gauss)

View File

@ -1,4 +1,3 @@
from collections import namedtuple
import logging
import math
import sys
@ -10,7 +9,7 @@ from hypothesis import given, assume, settings, HealthCheck
from hypothesis.strategies import floats, integers
import rasterio
from rasterio.errors import RasterioDeprecationWarning, WindowError
from rasterio.errors import WindowError
from rasterio.windows import (
crop, from_bounds, bounds, transform, evaluate, window_index, shape,
Window, intersect, intersection, get_data_window, union,
@ -631,3 +630,27 @@ def test_zero_height(sy):
"""Permit a zero height window"""
transform = Affine.translation(0, 45.0) * Affine.scale(1.0, sy)
assert from_bounds(0.0, 44.0, 1.0, 44.0, transform).height == 0
def test_union_boundless_left():
"""Windows entirely to the left of a dataset form a proper union"""
uw = union(
Window(col_off=-10, row_off=0, width=2, height=2),
Window(col_off=-8.5, row_off=0, width=2.5, height=2),
)
assert uw.col_off == -10
assert uw.width == 4
assert uw.height == 2
assert uw.row_off == 0
def test_union_boundless_above():
"""Windows entirely above a dataset form a proper union"""
uw = union(
Window(col_off=0, row_off=-10, width=2, height=2),
Window(col_off=0, row_off=-8.5, width=2, height=2.5),
)
assert uw.row_off == -10
assert uw.height == 4
assert uw.width == 2
assert uw.col_off == 0

View File

@ -107,8 +107,12 @@ def test_write_sbyte(tmpdir):
with rasterio.open(
name, 'w',
driver='GTiff', width=100, height=100, count=1,
dtype=a.dtype) as s:
s.write(a, indexes=1)
dtype=a.dtype) as dst:
dst.write(a, indexes=1)
with rasterio.open(name) as dst:
assert (dst.read() == -33).all()
info = subprocess.check_output(["gdalinfo", "-stats", name]).decode('utf-8')
assert "Minimum=-33.000, Maximum=-33.000, Mean=-33.000, StdDev=0.000" in info
assert 'SIGNEDBYTE' in info