RasterWizard

New in version 1.3.

class helpers.RasterWizard(layer=None)[source]

Get QGIS raster layers as numpy arrays and the result back into QGIS as a new raster layer.

To be used in the QGIS Python console for raster processing with numpy, SciPy, scikit-image, scikit-learn or other python libraries. Great for prototype development and experimenting with algorithms.

The resulting numpy array can be loaded back into QGIS as a new raster layer, as long as the number of pixels is the same as the input layer and the geotransform is not changed (no reprojection, no subsetting in numpy). The number of bands and the datatype can be different. See toarray() and tolayer().

On very large rasters, the processing can be done in windows (tiles) to avoid crashes, see get_windows().

Parameters

layer (QgsRasterLayer, optional) – instance of qgis.core.QgsRasterLayer, the layer to be processed. Optional, default is the active layer.

Example:

from scipy_filters.helpers import RasterWizard
from scipy import ndimage

wizard = RasterWizard() # Uses active layer if layer is not given
a = wizard.toarray()    # Returns numpy array with all bands

# Any calculation, for example a sobel filter with Scipy
# In the example, the result is a numpy array with dtype float32
b = ndimage.sobel(a, output="float32") 

# Write the result to a geotiff and load it back into QGIS
wizard.tolayer(b, name="Sobel", filename="/path/to/sobel.tif")

# You can also get pixel values at [x, y]
wizard[0,10] 

# or information like shape, CRS, etc.
wizard.shape   # like numpy_array.shape
wizard.crs_wkt # CRS as WKT string
wizard.crs     # CRS as QgsCoordinateReferenceSystem
property countbands

Number of bands in the raster layer.

Returns

Number of bands

Return type

int

property crs

Get the CRS of the layer as QgsCoordinateReferenceSystem.

Returns

CRS as qgis.core.QgsCoordinateReferenceSystem

Return type

QgsCoordinateReferenceSystem

property crs_id

Get the authority identifier of the CRS.

Returns a string like “EPSG:4326”.

Returns

CRS authority identifier

Return type

str

property crs_wkt

Get the CRS of the layer as WKT string.

Returns

CRS as WKT string

Return type

str

datatype(band=1, as_int=False)[source]

Get the name of the GDAL datatype of a band.

Parameters
  • band (int, optional) – Band index, default is 1. The fist band is indexed with 1.

  • as_int (bool, optional) – Return the datatype as integer instead, as defined in enum GDALDataType; default is False

Returns

GDAL Datatype as string or integer

property ds

GDAL dataset of the raster layer.

Returns

GDAL dataset

Return type

osgeo.gdal.Dataset

property filename

Path and filename of the raster layer source file.

Returns

Path and filename

Return type

str

property geotransform

Geotransform of the raster layer as used by GDAL.

The geotransform is a tuple with 6 elements of type float: top left x, w-e pixel resolution, rotation, top left y, rotation, n-s pixel resolution.

Returns

Geotransform

Return type

tuple

get_windows(windowsize=5000, margin=0)[source]

Generator to get windows for processing large rasters in tiles.

The generated windows are instances of RasterWindow. These do not contain the data, only the pixel indices and sizes that are used internally to read and write the data with GDAL. You can iterate over the generated windows, read the data with toarray() and write the data with write_window(). The output dataset should be set first with set_out_ds(), otherwise a virtual in-memory file is used. After processing all windows, the file is written and loaded back into QGIS with load_output().

Note that numpy, scipy etc. are very performant on large arrays. It is best to use a large windowsize or even the whole raster, as long as enough memory is avaible. The windows can have a margin, for algorthims that consider the neighborhood of a pixel as well. For example, a 3x3 kernel needs a margin of 1.

If you need the number of windows (e.g. for a progress bar), use number_of_windows().

Parameters
  • windowsize (int, optional) – Size of the windows in x and y direction in pixels, default is 5000

  • margin (int, optional) – Size of the margin in pixels, default is 0

Returns

Generator yielding RasterWindow instances

Example:

from scipy_filters.helpers import RasterWizard
from scipy import ndimage

wizard = RasterWizard()

wizard.set_out_ds(filename="/path/to/sobel.tif", dtype="float32")

for win in wizard.get_windows(windowsize=2048, margin=1):
    a = wizard.toarray(win=win)
    # Your calculation with numpy array a
    a = ndimage.sobel(a, output="float32")
    wizard.write_window(a, win)

wizard.load_output()
property has_axis_inverted

Returns True if the axis order of the CRS is inverted, i.e. not east/north (longitude/latitude).

Returns

True if axis order is inverted

Return type

bool

property is_geographic

Returns True if the CRS is geographic (i.e. not projected).

Returns

True if geographic CRS

Return type

bool

load_output(name='Wizard', stats=True)[source]

After processing windows, write the output file and load it back into QGIS.

The output file is written to the path given with set_out_ds() or to a virtual in-memory file. See get_windows() for a full example.

Parameters
  • name (str, optional) – Name of the new layer in QGIS, default is “Wizard”

  • stats (bool, optional) – Calculate and write band statistics (min, max, mean, std), default is True

Returns

QgsRasterLayer instance of the new layer

Return type

QgsRasterLayer

property map_units

Return the units for the projection used by the CRS.

Returns

Map units

Return type

Qgis.DistanceUnit enum

property name

Name of the raster layer in QGIS.

Returns

Layer name

Return type

str

property ndim

Number of dimensions of the raster layer.

Returns 2 for single band and 3 for multiple bands.

Returns

Number of dimensions

Return type

int

nodata(band=1)[source]

Get the no data value of a band. If no no data value is set, None is returned.

Note: In GeoTiff, the same no data value is used for all bands.

Parameters

band (int, optional) – Band index, default is 1. The fist band is indexed with 1.

Returns

No data value

Return type

float or None

number_of_windows(windowsize=5000)[source]

Returns the number of windows that would be generated with get_windows() with a given windowsize.

Since get_windows() is a generator, it can only be used once. Use this function to get the number of windows if you need it for a progress bar etc.

Parameters

windowsize (int, optional) – Size of the windows in pixels, default is 2048

Returns

Number of windows

Return type

int

property out_filename

Path and filename of the output geotiff file.

Is None if no output file is set. Normally, the output file is set when calling tolayer(). When iterating over windows, the output file should be set with set_out_ds().

Returns

Path or None

Return type

str or None

property pixel_height

Returns the height of a pixel in map units, read from the geotransform.

Returns

Pixel height

Return type

float

property pixel_width

Returns the width of a pixel in map units, read from the geotransform.

Returns

Pixel width

Return type

float

set_out_ds(filename=None, bands=None, dtype=None, nodata=None)[source]

Set up a new output dataset for writing.

Only to be used directly when iterating over windows, see get_windows() for example. Otherwise, the output is set in :py:meth`tolayer`.

Parameters
  • filename (str, optional) – Filename or full file path for the output geotiff, default is None (in-memory only)

  • bands (int, optional) – Number of bands in the output dataset, default is None (same as input)

  • dtype (str, optional) – Datatype of the output dataset, accepts numpy dtypes: “uint8”, “int8”, “uint16”, “int16”, “uint32”, “int32”, “float32”, “float64”, “complex64”, “complex128” default is None, using the datatype of the input layer (same as “input”).

  • nodata (int, optional) – No data value, default is None

property shape: tuple[int]

Shape of the raster layer

Given in the order (bands, x, y); identical to numpy_array.shape

Returns

Shape of the raster layer

Return type

tuple

toarray(band=None, win=None, wrapping=False, bands_last=False)[source]

Get the raster data as numpy array.

With parameter band, a 2D array is returned; otherwise, all bands are returned as a 3D array (or as 2D array if there is only one band). The order of the dimensions is [bands,] x, y.

Alternatively, by setting bands_last=True, the order of the dimensions is x, y [,bands] (e.g. expected by scikit-image).

Can be used together with the RasterWindow class to calculate in a moving window, see get_windows() for an example.

Parameters
  • band (int, optional) – Band index, default is None (all bands). The fist band is indexed with 1.

  • win (RasterWindow, optional) – RasterWindow instance for processing in a moving window, default is None

  • wrapping (bool, optional) – Fill the margins of the window with the data of the far side of the dataset, to be used with scipy filters with mode=’wrap’, default is False

  • bands_last (bool, optional) – If True, the order of the dimensions is x, y [,bands], default is False, expecting [bands,] x, y

Returns

Numpy array with the raster data. 2D if only one band, 3D if multiple bands.

Return type

numpy.ndarray

tolayer(array, name='Wizard', dtype='auto', filename=None, stats=True, bands_last=False, nodata=None)[source]

Load a numpy array as a new raster layer in QGIS.

The array must have the same dimensions in x and y as the input layer. The number of bands and the datatype can be different. First, the data is saved as geotiff if a filename (full file path) is given. Otherwise, the geotiff is in-memory only, using GDAL virtual file system.

If iterating over windows, use write_window() and load_output() instead.

Parameters
  • array (numpy.ndarray) – Numpy array with the raster data

  • name (str, optional) – Name of the new layer in QGIS, default is “Wizard”

  • dtype (str, optional) –

    Datatype of the output layer, accepts numpy dtypes: “uint8”, “int8”, “uint16”, “int16”, “uint32”, “int32”, “float32”, “float64”, “complex64”, “complex128” default is “auto”, using the GDAL datatype matching the numpy array datatype. “input” uses the datatype of the input layer.

  • filename (str, optional) – Filename or full file path for the output geotiff, default is None (in-memory only)

  • stats (bool, optional) – Calculate and write band statistics (min, max, mean, std), default is True

  • bands_last (bool, optional) – If True, the order of the dimensions is x, y [,bands], default is False, expecting [bands,] x, y

  • nodata (int, optional) – No data value, default is None

Returns

QgsRasterLayer instance of the new layer

Return type

QgsRasterLayer

write_window(array, win, band=None, bands_last=False)[source]

Write a numpy array into the window of the output dataset.

The window is an instance of RasterWindow, generated with get_windows(). The margin of the window is automatically sliced off before writing the data back to the raster. After processing all windows, the file is written and loaded back into QGIS with load_output(). For a full example with reading and writing windows and loading the result, see get_windows().

Parameters
  • array (numpy.ndarray) – Numpy array with the data to be written

  • win (RasterWindow) – Window, RasterWindow instance

  • band (int, optional) – Band index, default is None (all bands). The fist band is indexed with 1.

  • bands_last (bool, optional) – If True, the order of the dimensions is x, y [,bands], default is False, expecting [bands,] x, y