geofunctions

Objects

ProcessOptions

@dataclass
class ProcessOptions()

Data class to store all the options for processing.

Attributes:

append bool - append data to existing records by column, otherwise create new records. Defaults to True.
buffer float | int - buffer around centre to search pixels,
bounds List[float], optional - bbox of raster data will convert buffer value if input crs is 4326. Defaults to 1.
crs int, optional - input CRS to use for transforming.
date str - input date in YYYY-mm-DD format.
geotransform tuple - with transform matrix.
id_column str, optional - column that is used to identify unique geometries.
parameters str - name of the attribute describing a distance. Can be used instead of buffer to have a more dynamic search.
projection str - projection string in proj4 format.
simplify float | int - tolerance of simplification of geometry. Defaults to 0.
threshold float | int, optional - set a threshold which will cause the creation of a mask layer to mask values underneath the threshold value.
overwrite bool - overwrite existing records. Defaults to False.
XSize int - X dimension of tile in pixels.
YSize int - Y dimension of tile in pixels.
ALL_TOUCHED str - special GDAL parameter to regard all pixels touching the geometry. Defaults to True

get_bounds

def get_bounds(src: gdal.Dataset) -> List[float]

GIS function to calculate outer bounds of a raster dataset

Arguments:

src gdal.Dataset - Input raster dataset

Returns:

List[float] - bbox describing outer edge of raster dataset

add_raster_options

def add_raster_options(tile: gdal.Dataset | str)

Add specific GDAL raster options to ProcessOptions

Arguments:

tile gdal.Dataset | str - Input dataset or str pointer to raster file.

Functions

clean_polygons

def clean_polygons(polygons: List[dict], columns: dict | list) -> list

Filters data based on a whether records exists in the columns lookup table.

Arguments:

polygons List[dict] - input list of polygon records
columns dict | list - snippet of lookup table. Can also be an array of columns.

Returns:

list - filtered polygons

dissolve

def dissolve(geometries: List[str]) -> List[str]

GIS function to dissolve overlapping polygons into larger polygons. Attribute table is not guarded and one polygon will take precidence.

Arguments:

geometries List[str] - list of wkt strings

Returns:

List[str] - list of dissolved wkt strings length of list might have changed from input

Examples:

dissolve(
    geometries=[
        "POLYGON((804562.4672512376 688315.5459430797,804562.2603874448 688316.5859175948,804562.4150336159 688316.07611749,804562.4672512376 688315.5459430797))",
        "POLYGON((804535.3859031041 688334.3479520245,804535.1279857255 688335.6445902472,804535.3207982739 688335.0089724582,804535.3859031041 688334.3479520245))"
    ]
)

Notes:

Example geometries aren't necessarilly valid.

polygonize_array

def polygonize_array(data: np.ndarray,
                     transformation: tuple,
                     projection: str,
                     dn_value: int = 1) -> List[str]

GIS function to polygonize a numpy array to geometries

Arguments:

data np.ndarray - input data array.
transformation tuple - transformation matrix of feature.
projection str - input GDAL projection string src_ds.GetProjection().
dn_value int, optional - input burn value. Defaults to 1.

Returns:

List[str] - list of polygonized WKT strings

Examples:

polygonize_array(
    data=numpy.array([[1, 2, 3], [4, 5, 6]]),
    transformation=(399960.0, 30.0, 0.0, -2400000.0, 0.0, -30.0),
    projection='PROJCS["WGS 84 / UTM zone 21N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",-57],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32621"]]',
    dn_value=1
)

generate_bounding_box

def generate_bounding_box(envelope: tuple,
                          options: ProcessOptions) -> dict | None

Function to generate

Arguments:

envelope tuple - source envelope coordinates.
options ProcessOptions - processing options.

Returns:

dict| None: an object holding the transformation matrix and bbox

Examples:

generate_bounding_box(
    envelope=(1038796.0904789402, 1038798.0904789402, 1122432.7277924367, 1122434.7277924367),
    options=ProcessOptions(
        geotransform=(1038604.4586980281, 0.46867100000000056, 0.0, 1123215.2035654117, 0.0, -0.46867100000000894)
        XSize=4065,
        YSize=5000
    )
)

get_dtype

def get_dtype(polygon: dict) -> int

Estimate the numpy array dtype from GDAL or string representing a digit.

Arguments:

polygon dict - a polygon record.

Raises:

TypeError - raised if you try to burn funny types.

Returns:

int - the numpy dtype number.

zonal_stats

def zonal_stats(tile: gdal.Dataset,
                polygon: dict,
                options: ProcessOptions,
                attribute_name: str = "id") -> dict

GIS function to calculate zonal statistics from a PostGIS vector table

Arguments:

tile gdal.Dataset - input raster dataset.
polygon dict - input polygon information, representation is a row from a PostGIS table.
options ProcessOptions - input processing options.
attribute_name str, optional - primary key of vector table, or any other unique ID column. Defaults to "id".

Returns:

dict - output geometry object with polygonized WKT string and corresponding zonal stats data. Is empty if no significant pixels values were found.

Examples:

zonal_stats(
    tile=gdal.Dataset,
    polygon={
        "id": 1,
        "polygon": "POLYGON((804562.4672512376 688315.5459430797,804562.2603874448 688316.5859175948,804562.4150336159 688316.07611749,804562.4672512376 688315.5459430797))",
        "date": "2021-01-01",
        "meandist": 0.5
    },
    geotransform=(399960.0, 30.0, 0.0, -2400000.0, 0.0, -30.0),
    projection='PROJCS["WGS 84 / UTM zone 21N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",-57],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32621"]]',
    options=ProcessOptions(foo="bar"),
    attribute_name=id_0
)

zonal_stats_polygons

def zonal_stats_polygons(tile: gdal.Dataset | str, polygons: List[dict],
                         rasterlayer: str,
                         options: ProcessOptions) -> List[dict]

Collector function that cleanses the incoming list of polygons and runs the zonal statistics per polygon. It also parses the zonal statistics return into a uncluttered array of objects. Function also checks which geometries might have been missed.

Arguments:

tile gdal.Dataset | str - Input dataset or str pointer to raster file.
polygons List[dict] - Input list of polygon objects as returned from a PostGIS table.
rasterlayer str - Input name of source raster layer.
options ProcessOptions - Input processing options.

Raises:

Exception - panics on generic exceptions.

Returns:

List[dict] - array of polygon objects

Examples:

zonal_stats_polygons(
    tile="/mnt/data/some.tiff",
    polygons=
    [
        {
            "id": 1,
            "polygon": ("POLYGON(
                (
                    804562.4672512376 688315.5459430797,
                    804562.4150336159 688316.07611749,
                    804562.4167345229 688316.076381122,
                    804562.4672512376 688315.5459430797
                )
            )"),
            "date": "2021-01-01",
            "meandist": 0.5
        }
    ],
    rasterlayer="objectheightmodel",
    options={
        "foo": "bar"
    }
)

generate_bin_mask

def generate_bin_mask(
    tile: gdal.Dataset, threshold: float | int,
    bbox: list | tuple = ()) -> np.ndarray

Function to generate a binary mask based on a threshold value

Arguments:

tile gdal.Dataset - input raster data.
threshold float | int - the threshold value.
bbox list | tuple, optional - an optional bbox, if not given it will use the entire array. Defaults to ().

Returns:

np.ndarray - masked array of the raster data.

Objects​

ProcessOptions​

get_bounds​

add_raster_options​

Functions​

clean_polygons​

dissolve​

polygonize_array​

generate_bounding_box​

get_dtype​

zonal_stats​

zonal_stats_polygons​

generate_bin_mask​

Objects

ProcessOptions

get_bounds

add_raster_options

Functions

clean_polygons

dissolve

polygonize_array

generate_bounding_box

get_dtype

zonal_stats

zonal_stats_polygons

generate_bin_mask