geoutils.georaster.raster

geoutils.georaster provides a toolset for working with raster data.

Contents

Classes

Raster

class geoutils.georaster.raster.Raster(filename_or_dataset, bands=None, load_data=True, downsample=1, masked=True, nodata=None, attrs=None, as_memfile=False)[source]

Bases: object

Create a Raster object from a rasterio-supported raster dataset.

If not otherwise specified below, attribute types and contents correspond to the attributes defined by rasterio.

Attributes:
filenamestr

The path/filename of the loaded, file, only set if a disk-based file is read in.

datanp.array

Loaded image. Dimensions correspond to (bands, height, width).

nbandsint

Number of bands loaded into .data

bandstuple

The indexes of the opened dataset which correspond to the bands loaded into data.

is_loadedbool

True if the image data have been loaded into this Raster.

dsrio.io.DatasetReader

Link to underlying DatasetReader object.

bounds

count

crs

dataset_mask

driver

dtypes

height

indexes

name

nodata

res

shape

transform

width

__init__(filename_or_dataset, bands=None, load_data=True, downsample=1, masked=True, nodata=None, attrs=None, as_memfile=False)[source]

Load a rasterio-supported dataset, given a filename.

Parameters

  • filename_or_dataset – The filename of the dataset.

  • bands – The band(s) to load into the object. Default is to load all bands.

  • load_data – Load the raster data into the object. Default is True.

  • downsample – Reduce the size of the image loaded by this factor. Default is 1

  • masked – the data is loaded as a masked array, with no data values masked. Default is True.

  • nodata – nodata to be used (overwrites the metadata). Default is None, i.e. reads from metadata.

  • attrs – Additional attributes from rasterio’s DataReader class to add to the Raster object. Default list is [‘bounds’, ‘count’, ‘crs’, ‘dataset_mask’, ‘driver’, ‘dtypes’, ‘height’, ‘indexes’, ‘name’, ‘nodata’, ‘res’, ‘shape’, ‘transform’, ‘width’] - if no attrs are specified, these will be added.

  • as_memfile – open the dataset via a rio.MemoryFile.

Returns

A Raster object

astype(dtype: np.dtype | type | str, inplace: typing_extensions.Literal[False]) → geoutils.georaster.raster.Raster[source]
astype(dtype: np.dtype | type | str, inplace: typing_extensions.Literal[True]) → None

Converts the data type of a Raster object.

Parameters

  • dtype – Any numpy dtype or string accepted by numpy.astype

  • inplace – Set to True to modify the raster in place.

Returns

the output Raster with dtype changed.

coords(offset='corner', grid=True)[source]

Get x,y coordinates of all pixels in the raster.

Parameters

  • offset – coordinate type. If ‘corner’, returns corner coordinates of pixels. If ‘center’, returns center coordinates. Default is corner.

  • grid – Return grid

Returns x,y

numpy arrays corresponding to the x,y coordinates of each pixel.

copy(new_array=None)[source]

Copy the Raster object in memory

Parameters

new_array – New array to use for the copied Raster

Returns

crop(cropGeom, mode='match_pixel', inplace=True)[source]

Crop the Raster to a given extent.

Parameters

  • cropGeom – Geometry to crop raster to, as either a Raster object, a Vector object, or a list of coordinates. If cropGeom is a Raster, crop() will crop to the boundary of the raster as returned by Raster.ds.bounds. If cropGeom is a Vector, crop() will crop to the bounding geometry. If cropGeom is a list of coordinates, the order is assumed to be [xmin, ymin, xmax, ymax].

  • mode – one of ‘match_pixel’ (default) or ‘match_extent’. ‘match_pixel’ will preserve the original pixel resolution, cropping to the extent that most closely aligns with the current coordinates. ‘match_extent’ will match the extent exactly, adjusting the pixel resolution to fit the extent.

  • inplace – Update the raster inplace or return copy.

Returns

None if inplace=True and a new Raster if inplace=False

property data: np.ndarray | np.ma.masked_array

Get data.

Returns

data array.

classmethod from_array(data, transform, crs, nodata=None)[source]

Create a Raster from a numpy array and some geo-referencing information.

Parameters

  • data – data array

  • transform – the 2-D affine transform for the image mapping. Either a tuple(x_res, 0.0, top_left_x, 0.0, y_res, top_left_y) or an affine.Affine object.

  • crs – Coordinate Reference System for image. Either a rasterio CRS, or the EPSG integer.

  • nodata – nodata value

Returns

A Raster object containing the provided data.

Example:

You have a data array in EPSG:32645. It has a spatial resolution of 30 m in x and y, and its top left corner is X=478000, Y=3108140.

>>> transform = (30.0, 0.0, 478000.0, 0.0, -30.0, 3108140.0)
>>> myim = Raster.from_array(data, transform, 32645)
get_bounds_projected(out_crs, densify_pts_max=5000)[source]

Return self’s bounds in the given CRS.

Parameters

  • out_crs (CRS) – Output CRS

  • densify_pts_max (int) – Maximum points to be added between image corners to account for non linear edges. Reduce if time computation is really critical (ms) or increase if extent is not accurate enough.

Return type

BoundingBox

ij2xy(i, j, offset='center')[source]

Return x,y coordinates for a given row, column index pair.

Parameters

  • i – row (i) index of pixel.

  • j – column (j) index of pixel.

  • offset – return coordinates as “corner” or “center” of pixel

Returns x, y

x,y coordinates of i,j in reference system.

info(stats=False)[source]

Returns string of information about the raster (filename, coordinate system, number of columns/rows, etc.).

Parameters

stats (bool) – Add statistics for each band of the dataset (max, min, median, mean, std. dev.). Default is to not calculate statistics.

Return type

str

Returns

text information about Raster attributes.

interp_points(pts, input_latlon=False, mode='linear', band=1, area_or_point=None, **kwargs)[source]

Interpolate raster values at a given point, or sets of points.

Parameters

pts – Point(s) at which to interpolate raster value. If points fall outside of image,

value returned is nan. Shape should be (N,2)’ :param input_latlon: Whether the input is in latlon, unregarding of Raster CRS :param mode: One of ‘linear’, ‘cubic’, or ‘quintic’. Determines what type of spline is

used to interpolate the raster value at each point. For more information, see scipy.interpolate.interp2d. Default is linear.

Parameters

  • band – Raster band to use

  • area_or_point – shift index according to GDAL AREA_OR_POINT attribute (None) or force position (‘Point’ or ‘Area’) of the interpretation of where the raster value corresponds to in the pixel (‘Area’ = lower left or ‘Point’ = center)

Returns rpts

Array of raster value(s) for the given points.

intersection(rst)[source]

Returns the bounding box of intersection between this image and another.

If the rasters have different projections, the intersection extent is given in self’s projection system.

:param rst : path to the second image (or another Raster instance)

Returns

extent of the intersection between the 2 images (xmin, ymin, xmax, ymax) in self’s coordinate system.

property is_modified: bool

Check whether file has been modified since it was created/opened.

Return type

bool

Returns

True if Raster has been modified.

load(bands=None, **kwargs)[source]

Load specific bands of the dataset, using rasterio.read(). Ensure that self.data.ndim = 3 for ease of use (needed e.g. in show)

Parameters

bands – The band(s) to load. Note that rasterio begins counting at 1, not 0.

**kwargs: any additional arguments to rasterio.io.DatasetReader.read. Useful ones are: .. hlist:: * out_shape : to load a subsampled version * window : to load a cropped version * resampling : to set the resampling algorithm

outside_image(xi, yj, index=True)[source]

Check whether a given point falls outside of the raster.

Parameters

  • xi (ndarray) – Indices (or coordinates) of x direction to check.

  • yj (ndarray) – Indices (or coordinates) of y direction to check.

  • index (bool) – Interpret ij as raster indices (default is True). If False, assumes ij is coordinates.

Returns is_outside

True if ij is outside of the image.

Return type

bool

polygonize(in_value=1)[source]

Return a GeoDataFrame polygonized from a raster.

Parameters

in_value – Value or range of values of the raster from which to create geometries (Default is 1). If ‘all’, all unique pixel values of the raster are used.

Returns

Vector containing the polygonized geometries.

reproject(dst_ref=None, dst_crs=None, dst_size=None, dst_bounds=None, dst_res=None, dst_nodata=None, src_nodata=None, dtype=None, resampling=Resampling.nearest, silent=False, n_threads=0, memory_limit=64)[source]

Reproject raster to a specified grid.

The output grid can either be given by a reference Raster (using dst_ref), or by manually providing the output CRS (dst_crs), dimensions (dst_size), resolution (with dst_size) and/or bounds (dst_bounds). Any resampling algorithm implemented in rasterio can be used.

Currently: requires image data to have been loaded into memory. NOT SUITABLE for large datasets yet! This requires work…

To reproject a Raster with different source bounds, first run Raster.crop.

Parameters

  • dst_ref – a reference raster. If set will use the attributes of this raster for the output grid. Can be provided as Raster/rasterio data set or as path to the file.

  • crs – Specify the Coordinate Reference System to reproject to. If dst_ref not set, defaults to self.crs.

  • dst_size – Raster size to write to (x, y). Do not use with dst_res.

  • dst_bounds – a BoundingBox object or a dictionary containing left, bottom, right, top bounds in the source CRS.

  • dst_res – Pixel size in units of target CRS. Either 1 value or (xres, yres). Do not use with dst_size.

  • dst_nodata – nodata value of the destination. If set to None, will use the same as source, and if source is None, will use GDAL’s default.

  • src_nodata – nodata value of the source. If set to None, will read from the metadata.

  • resampling – A rasterio Resampling method

  • silent – If True, will not print warning statements

  • n_threads – The number of worker threads. Defaults to (os.cpu_count() - 1).

  • memory_limit – The warp operation memory limit in MB. Larger values may perform better.

Returns

Raster

save(filename, driver='GTiff', dtype=None, compress='deflate', tiled=False, blank_value=None, co_opts=None, metadata=None, gcps=None, gcps_crs=None)[source]

Write the Raster to a geo-referenced file.

Given a filename to save the Raster to, create a geo-referenced file on disk which contains the contents of self.data.

If blank_value is set to an integer or float, then instead of writing the contents of self.data to disk, write this provided value to every pixel instead.

Parameters

  • filename – Filename to write the file to.

  • driver – the ‘GDAL’ driver to use to write the file as.

  • dtype – Data Type to write the image as (defaults to dtype of image data)

  • compress – Compression type. Defaults to ‘deflate’ (equal to GDALs: COMPRESS=DEFLATE)

  • tiled – Whether to write blocks in tiles instead of strips. Improves read performance on large files, but increases file size.

  • blank_value – Use to write an image out with every pixel’s value corresponding to this value, instead of writing the image data to disk.

  • co_opts – GDAL creation options provided as a dictionary, e.g. {‘TILED’:’YES’, ‘COMPRESS’:’LZW’}

  • metadata – pairs of metadata key, value

  • gcps – list of gcps, each gcp being [row, col, x, y, (z)]

  • gcps_crs – the CRS of the GCPS (Default is None)

Returns

None.

set_ndv(ndv, update_array=False)[source]

Set new nodata values for bands (and possibly update arrays).

Parameters

  • ndv – nodata values

  • update_array – change the existing nodata in array

shift(xoff, yoff)[source]

Translate the Raster by a given x,y offset.

Parameters

  • xoff (float) – Translation x offset.

  • yoff (float) – Translation y offset.

Return type

None

show(band=None, cmap=None, vmin=None, vmax=None, cb_title=None, add_cb=True, ax=None, **kwargs)[source]

Show/display the image, with axes in projection of image.

This method is a wrapper to rasterio.plot.show. Any **kwargs which you give this method will be passed to rasterio.plot.show.

Parameters

  • band – which band to plot, from 0 to self.count-1 (default is all)

  • cmap – The figure’s colormap. Default is plt.rcParams[‘image.cmap’]

  • vmin – Colorbar minimum value. Default is data min.

  • vmax – Colorbar maximum value. Default is data min.

  • cb_title – Colorbar label. Default is None.

  • add_cb – Set to True to display a colorbar. Default is True.

  • ax – A figure ax to be used for plotting. If None, will create default figure and axes,and plot figure directly.

Returns

if ax is not None, returns (ax, cbar) where cbar is the colorbar (None if add_cb is False)

You can also pass in **kwargs to be used by the underlying imshow or contour methods of matplotlib. The example below shows provision of a kwarg for rasterio.plot.show, and a kwarg for matplotlib as well:

import matplotlib.pyplot as plt
ax1 = plt.subplot(111)
mpl_kws = {'cmap':'seismic'}
myimage.show(ax=ax1, mpl_kws)
split_bands(copy=False, subset=None)[source]

Split the bands into separate copied rasters.

Parameters

  • copy – Copy the bands or return slices of the original data.

  • subset – Optional. A subset of band indices to extract. Defaults to all.

Returns

A list of Rasters for each band.

to_points(subset: float | int, as_frame: typing_extensions.Literal[True], pixel_offset: typing_extensions.Literal[center, corner]) → geopandas.geodataframe.GeoDataFrame[source]
to_points(subset: float | int, as_frame: typing_extensions.Literal[False], pixel_offset: typing_extensions.Literal[center, corner]) → numpy.ndarray

Subset a point cloud of the raster.

If ‘subset’ is either 1 or is equal to the pixel count, all points are returned in order. If ‘subset’ is smaller than 1 (for fractions) or the pixel count, a random sample is returned.

If the raster is not loaded, sampling will be done from disk without loading the entire Raster.

Formats:

  • as_frame == None | False: A numpy ndarray of shape (N, 2 + nbands) with the columns [x, y, b1, b2..].

  • as_frame == True: A GeoPandas GeoDataFrame with the columns [“b1”, “b2”, …, “geometry”]

Parameters

  • subset – The point count or fraction. If ‘subset’ > 1, it’s parsed as a count.

  • as_frame – Return a GeoDataFrame with a geometry column and crs instead of an ndarray.

  • pixel_offset – The point at which to associate the pixel coordinate with (‘corner’ == upper left).

Raises

ValueError – If the subset count or fraction is poorly formatted.

Returns

An ndarray/GeoDataFrame of the shape (N, 2 + nbands) where N is the subset count.

to_xarray(name=None)[source]

Convert this Raster into an xarray DataArray using rioxarray.

This method uses rioxarray to generate a DataArray with associated geo-referencing information.

See the documentation of rioxarray and xarray for more information on the methods and attributes of the resulting DataArray.

Parameters

name – Set the name of the DataArray.

Returns

xarray DataArray

value_at_coords(x, y, latlon=False, band=None, masked=False, window=None, return_window=False, boundless=True, reducer_function=<numpy.ma.core._frommethod object>)[source]

Extract the pixel value(s) at the nearest pixel(s) from the specified coordinates.

Extract pixel value of each band in dataset at the specified coordinates. Alternatively, if band is specified, return only that band’s pixel value.

Optionally, return mean of pixels within a square window.

Parameters

  • x – x (or longitude) coordinate.

  • y – y (or latitude) coordinate.

  • latlon – Set to True if coordinates provided as longitude/latitude.

  • band – the band number to extract from.

  • masked – If masked is True the return value will be a masked array. Otherwise (the default) the return value will be a regular array.

  • window – expand area around coordinate to dimensions window * window. window must be odd.

  • return_window – If True when window=int, returns (mean,array) where array is the dataset extracted via the specified window size.

  • boundless – If True, windows that extend beyond the dataset’s extent are permitted and partially or completely filled arrays (with self.nodata) will be returned as appropriate.

  • reducer_function – a function to apply to the values in window.

Returns

When called on a Raster or with a specific band set, return value of pixel.

Returns

If multiple band Raster and the band is not specified, a dictionary containing the value of the pixel in each band.

Returns

In addition, if return_window=True, return tuple of (values, arrays)

Examples

>>> self.value_at_coords(-48.125,67.8901,window=3)
Returns mean of a 3*3 window:
    v v v             v c v  | = float(mean)
    v v v /
(c = provided coordinate, v= value of surrounding coordinate)
xy2ij(x, y, op=<class 'numpy.float32'>, area_or_point=None, precision=None)[source]

Return row, column indices for a given x,y coordinate pair.

Parameters

  • x – x coordinates

  • y – y coordinates

  • op – operator to calculate index

  • precision – precision for rio.Dataset.index

  • area_or_point – shift index according to GDAL AREA_OR_POINT attribute (None) or force position (‘Point’ or ‘Area’) of the interpretation of where the raster value corresponds to in the pixel (‘Area’ = lower left or ‘Point’ = center)

Returns i, j

indices of x,y in the image.