Classifying Data Spatially

mapteksdk 2.0+

You can classify points based on their location relative to another entity, such as a surface A two-dimensional structure in 3D space composed of points connected by triangular facets or plane, using a spatial filter. For example, you can find the points in an object that lie entirely above a surface and then hide them.

The mapteksdk.data.spatial_filters subpackage provides functions to create various kinds of spatial filters. Each function returns a SpatialFilter object that can be applied to a set of points or combined with other b.

The main filters are summarised in the table below.

Spatial filter Matches points that lie...

above_surface()

above the specified surface

below_surface()

below the specified surface

on_surface()

exactly on the specified surface

inside_surface_extent()

within the surface’s boundary when viewed from above (ignores vertical position)

outside_surface_extent()

outside the surface’s boundary when viewed from above (ignores vertical position)

inside_solid()

within the specified solid

outside_solid()

outside the specified solid

above_plane()

above the specified plane

below_plane()

below the specified plane

on_plane()

exactly on the specified plane

inside_object_extent()

within an object’s bounding box The smallest axis-aligned, rectangular prism that completely contains the geometry of an object.

outside_object_extent()

outside an object’s bounding box The smallest axis-aligned, rectangular prism that completely contains the geometry of an object.

Below you will find examples of how to use various spatial filters. Most of these examples use the spatial_filter_demo.maptekobj dataset, which you can download to try out the scripts.

Classifying points using spatial filters

After creating a spatial filter, you can use it to classify points by calling the filter’s SpatialFilter.filter() method. The method takes the object or points to classify as an argument, and returns a Boolean array where each element is True if the corresponding point matches the filter, and False otherwise. For example, if you call filter() on the above_surface() filter and the first point in the supplied set lies above the surface, the first value of the array will be True. Once you have the Boolean array, you can use it to hide, colour, or otherwise process the points that match the filter.

Above or below surface

To classify points as above or below a surface, use the above_surface() and below_surface() spatial filter functions.

The following script prompts the user to pick An action in which the user selects an object, primitive, or coordinate, typically by clicking in a view window. a surface and an object to filter. It then hides all points (or blocks) in the object that lie above the picked surface, using the Boolean array returned by filter() as a mask on the point (or block) visibility array.

from mapteksdk.data import Surface, spatial_filters
from mapteksdk.operations import object_pick
from mapteksdk.project import Project

def main(project: Project):
    surface_oid = object_pick(
        object_types=(Surface), label="Pick surface to filter blocks or points above"
    )
    object_oid = object_pick(label="Pick object to filter")

    spatial_filter = spatial_filters.above_surface(surface_oid)

    with project.edit(object_oid) as data_object:
        matched_primitives = spatial_filter.filter(data_object)


        if hasattr(data_object, "point_visibility"):
            data_object.point_visibility[matched_primitives] = False
        elif hasattr(data_object, "block_visibility"):
            data_object.block_visibility[matched_primitives] = False
        else:
            # This should be unreachable, as the call to `filter()` should
            # have failed in this case.
            raise RuntimeError("Unable to filter object!")

if __name__ == "__main__":
    with Project() as main_project:
        with main_project.undo():
            main(main_project)

The following animations demonstrate running this script on a couple of different object types.

Running the script on a block model.

Running the script on a point set.

Note that in the block model example, some unfiltered blocks protrude above the topography after filtering is applied. This occurs because spatial filters operate on block centroids rather than the full block extent. As long as a block’s centroid lies below the surface, it is not filtered.

You can modify the script to hide points that lie below the surface by using below_surface() in place of above_surface(), as follows:

spatial_filter = spatial_filters.below_surface(surface_oid)

The following animation demonstrates running the script incorporating this change on the same point set:

Colouring matched points

Instead of hiding points matched by the filter, you can change their colour to make them stand out. The following line of code sets the colour of matched points to yellow by using the Boolean array as a mask on the point colour array:

data_object.point_colours[matched_primitives] = [200, 200, 0, 255]

The following animation demonstrates running the script with this line in place of the point-hiding code:

Storing classifications as a point attribute

You can also store the classifications as a point attribute on the object for later processing, as demonstrated in the following snippet:

# Points are not above the surface by default.
data_object.point_attributes["above"] = False
data_object.point_attributes["above"][matched_primitives] = True

This snippet first initialises all attribute values to False, then sets the attribute value to True only for those points that match the filter.

Including points that lie exactly on the surface

By default, points that lie exactly on a surface are considered to be neither above nor below it. To include such points, set include_on=True in above_surface() or below_surface().

The following snippet creates a filter that will match points that lie above or exactly on the surface:

spatial_filter = spatial_filters.above_surface(surface_oid, include_on=True)

Classifying points as above or below facets

If your script is working on an array or some subset of facets instead of a surface object, you can use the above_facets() and below_facets() functions to classify points relative to those facets.

The following snippet creates a filter that will match points that lie below a set of facets:

spatial_filter = spatial_filters.below_facets(points, facets)

Inside or outside a surface extent

To classify points as inside or outside a surface extent, use the inside_surface_extent() and outside_surface_extent() spatial filter functions. Points within the surface boundary when viewed from above are considered to be inside, while points outside the boundary are considered outside. In both cases, the point’s vertical position (above, below, or on the surface) is ignored.

The following snippet creates a filter that will match points that lie outside a surface boundary:

spatial_filter = spatial_filters.outside_surface_extent(surface_oid)

Inside or outside a solid

To classify points as inside or outside a solid A closed 3D surface that fully encloses a volume., use the inside_solid() and outside_solid() spatial filter functions.

Note:  Unlike the above and below surface filters, points that lie on a solid are considered to be inside the solid.

The following script is a simple modification of the script introduced above. The script prompts for a solid and a block model. It then colours yellow all blocks in the block model that lie within the solid.

from mapteksdk.data import Surface, spatial_filters
from mapteksdk.operations import object_pick
from mapteksdk.project import Project

def main(project: Project):
    surface_oid = object_pick(
        object_types=(Surface), label="Pick surface to filter blocks above"
    )
    object_oid = object_pick(label="Pick object to filter")

    spatial_filter = spatial_filters.inside_solid(surface_oid)

    with project.edit(object_oid) as data_object:
        matched_primitives = spatial_filter.filter(data_object)


        if hasattr(data_object, "point_colours"):
            data_object.point_colours[matched_primitives] = [220, 220, 0, 255]
        elif hasattr(data_object, "block_colours"):
            data_object.block_colours[matched_primitives] = [220, 220, 0, 255]
        else:
            # This should be unreachable, as the call to `filter()` should
            # have failed in this case.
            raise RuntimeError("Unable to filter object!")

if __name__ == "__main__":
    with Project() as main_project:
        with main_project.undo():
            main(main_project)

The following animation demonstrates running this script with a solid and block model from the example dataset. A section view is enabled to allow the solid, which lies wholly inside the block model, to be picked.

Note that some uncoloured blocks extend inside the solid, and some coloured blocks extend outside it. This occurs because the spatial filter operates on the block centroid rather than the full block extent—any block whose centre lies within the solid is considered a match.

The outside_solid() spatial filter can be used to perform the opposite operation—matching all points or blocks that lie outside the supplied solid. A simple change to the script above will cause the blocks outside the solid to be coloured instead of the blocks inside:

spatial_filter = spatial_filters.outside_solid(surface_oid)

The following animation demonstrates running the script with this change incorporated:

Note

Above or below a plane

To classify points as above or below a plane, use the above_plane() and below_plane() spatial filter functions. These are similar to the above_surface() and below_surface() filters, except they classify points with respect to a plane rather than a surface.

Note:  Whether a point is classified as above or below a plane is determined by the plane’s normal vector, not by its absolute height. Points on the side the normal points toward are considered above, while points on the opposite side are considered below. Inverting the normal vector reverses which side is treated as above.

The following script prompts the user to select an object and colours points above an XY plane passing through the object’s centre yellow.

from mapteksdk.data import Topology, spatial_filters
from mapteksdk.geometry import Plane
from mapteksdk.operations import object_pick
from mapteksdk.project import Project

def main(project: Project):
    object_oid = object_pick(object_types=(Topology,), label="Pick object to filter")

    with project.edit(object_oid) as data_object:
        spatial_filter = spatial_filters.above_plane(
            Plane.xy(z=data_object.extent.centre[2])
        )
        matched_primitives = spatial_filter.filter(data_object)

        if hasattr(data_object, "point_colours"):
            data_object.point_colours[matched_primitives] = [220, 220, 0, 255]
        elif hasattr(data_object, "block_colours"):
            data_object.block_colours[matched_primitives] = [220, 220, 0, 255]
        else:
            # This should be unreachable, as the call to `filter()` should
            # have failed in this case.
            raise RuntimeError("Unable to filter object!")

if __name__ == "__main__":
    with Project() as main_project:
        with main_project.undo():
            main(main_project)

The following animation demonstrates running this script on a block model:

Including points that lie exactly on the plane

By default, points that lie exactly on the plane are considered to be neither above nor below it. To include such points, set include_on=True in above_plane() and below_plane() .

For example, the following filter will match points that lie above or exactly on the XY plane through the centre of the object:

spatial_filter = spatial_filters.above_plane(
    Plane.xy(data_object.extent.centre[2]),
    include_on=True
)

Inside or outside an extent

Sometimes you only need to know whether points lie within or outside an object’s overall extent The set of minimum and maximum X, Y, and Z coordinates that define the axis-aligned bounding box enclosing an object., not its exact volume. The inside_object_extent() and outside_object_extent() filters match points within or outside an object’s axis-aligned bounding box—the smallest box that fully contains the object’s geometry. The following script demonstrates using this filter to colour an object.

from mapteksdk.data import Topology, spatial_filters
from mapteksdk.operations import object_pick
from mapteksdk.project import Project

def main(project: Project):
    extent_source = object_pick(
        object_types=(Topology), label="Pick object to provide the extent"
    )
    object_oid = object_pick(object_types=(Topology,), label="Pick object to filter")

    with project.edit(object_oid) as data_object:
        spatial_filter = spatial_filters.inside_object_extent(
            extent_source
        )
        matched_primitives = spatial_filter.filter(data_object)

        if hasattr(data_object, "point_colours"):
            data_object.point_colours[matched_primitives] = [220, 220, 0, 255]
        elif hasattr(data_object, "block_colours"):
            data_object.block_colours[matched_primitives] = [220, 220, 0, 255]
        else:
            # This should be unreachable, as the call to `filter()` should
            # have failed in this case.
            raise RuntimeError("Unable to filter object!")

if __name__ == "__main__":
    with Project() as main_project:
        with main_project.undo():
            main(main_project)

The following animation demonstrates running this script on a solid located within a block model. A section view is enabled to allow the solid, which lies wholly inside the block model, to be picked. The script colours the blocks within the solid’s bounding box yellow. The bounding box is an axis-aligned rectangular prism.

We can change the script to colour blocks outside the solid’s extents by using the outside_object_extent() filter in place of the inside_object_extent() function:

spatial_filter = spatial_filters.outside_object_extent(extent_source)

The following screenshot demonstratesthe result when the script is run with this change incorporated:

Using a user-defined extent

Instead of filtering by the extent of an object, you can also classify by a custom extent you define using the inside_extent() spatial filter. The following script demonstrates selecting an object and colouring points within a hard-coded extent:

from mapteksdk.data import Topology, spatial_filters
from mapteksdk.geometry import Extent
from mapteksdk.operations import object_pick
from mapteksdk.project import Project

def main(project: Project):
    object_oid = object_pick(object_types=(Topology,), label="Pick object to filter")

    with project.edit(object_oid) as data_object:
        spatial_filter = spatial_filters.inside_extent(
            Extent(
                (740400, 7642500, 150),
                (741000, 7643000, 200)
            )
        )
        matched_primitives = spatial_filter.filter(data_object)

        if hasattr(data_object, "point_colours"):
            data_object.point_colours[matched_primitives] = [220, 220, 0, 255]
        elif hasattr(data_object, "block_colours"):
            data_object.block_colours[matched_primitives] = [220, 220, 0, 255]
        else:
            # This should be unreachable, as the call to `filter()` should
            # have failed in this case.
            raise RuntimeError("Unable to filter object!")

if __name__ == "__main__":
    with Project() as main_project:
        with main_project.undo():
            main(main_project)

This matches this very small section of the block model, as demonstrated in the screenshot below.

To select points outside the same region instead, use the outside_extent() filter by replacing the line in the script that uses inside_extent() in the script with the following:

spatial_filter = spatial_filters.outside_extent(
    Extent(
        (740400, 7642500, 150),
        (741000, 7643000, 200)
    )
)

This matches the inverse of the blocks previously selected, as shown below.

Combining and inverting filters

Spatial filters can be combined or inverted to create more complex filtering conditions. You can use the & (AND) and | (OR) operators to combine filters, invert a filter to select points or blocks not matched by it, or use helper functions like all_filters() and any_filters() to work with multiple sub-filters programmatically.

Note:  Spatial filter objects are immutable once created. Combining or inverting filters does not change the originals; instead, new filters that represent the combined or inverted condition are returned.

Combining two filters with a logical AND

Two spatial filters can be combined with the & operator to create a new filter that matches only points matched by both. For example, the following script combines an above_surface() and a below_surface() filter to colour points located between a bottom and top surface picked by the user:

from mapteksdk.data import Surface, spatial_filters
from mapteksdk.operations import object_pick
from mapteksdk.project import Project

def main(project: Project):
    bottom_surface_oid = object_pick(
        object_types=(Surface), label="Pick the bottom surface"
    )
    top_surface_oid = object_pick(object_types=(Surface), label="Pick the top surface")
    object_oid = object_pick(label="Pick object to filter")

    spatial_filter = spatial_filters.above_surface(
        bottom_surface_oid
    ) & spatial_filters.below_surface(top_surface_oid)

    with project.edit(object_oid) as data_object:
        matched_primitives = spatial_filter.filter(data_object)

        if hasattr(data_object, "point_colours"):
            data_object.point_colours[matched_primitives] = [220, 220, 0, 255]
        elif hasattr(data_object, "block_colours"):
            data_object.block_colours[matched_primitives] = [220, 220, 0, 255]
        else:
            # This should be unreachable, as the call to `filter()` should
            # have failed in this case.
            raise RuntimeError("Unable to filter object!")

if __name__ == "__main__":
    with Project() as main_project:
        with main_project.undo():
            main(main_project)

The following animation shows the script in action. It uses a translated copy of a topography to provide the top surface:

All spatial filters support being combined with the & operator.

Combining multiple filters using a logical AND

The all_filters() function combines multiple spatial filters into a single filter that matches points or blocks only if all of the sub-filters match—essentially performing a logical AND across the group.

The following snippet demonstrates combining three spatial filters into a single filter using all_filters():

inside_all_filter = spatial_filters.all_filters(
    spatial_filters.inside_solid(surface_a),
    spatial_filters.inside_solid(surface_b),
    spatial_filters.inside_solid(surface_c),
)

The same result could be achieved by combining the filters using the & operator, as follows:

inside_all_filter = (
    spatial_filters.inside_solid(surface_a) &
    spatial_filters.inside_solid(surface_b) &
    spatial_filters.inside_solid(surface_c)
)

The main advantage of using all_filters() over the & operator is that it allows sub-filters to be created dynamically. For instance, the following snippet uses list comprehension to generate a list of inside_solid() filters, which are then passed to all_filters() via unpacking:

surface_ids: Sequence[ObjectID[Surface]]
inside_all_filter = spatial_filters.all_filters(
    *[
        spatial_filters.inside_solid(surface_id)
        for surface_id in surface_ids
    ]
)

This is a more convenient way of combining filters when the number of filters is not known in advance.

The animation below demonstrates the script in action. Only blocks located inside all three yellow cubes are coloured yellow; all other blocks are coloured cyan:

Combining two filters using a logical OR

Two spatial filters can be combined with the | operator to create a new filter that matches points matched by either filter. For example, the following script colours points yellow that are either above the top surface or below the bottom surface:

from mapteksdk.data import Surface, spatial_filters
from mapteksdk.operations import object_pick
from mapteksdk.project import Project

def main(project: Project):
    bottom_surface_oid = object_pick(
        object_types=(Surface), label="Pick the bottom surface"
    )
    top_surface_oid = object_pick(object_types=(Surface), label="Pick the top surface")
    object_oid = object_pick(label="Pick object to filter")

    spatial_filter = spatial_filters.above_surface(
        top_surface_oid
    ) | spatial_filters.below_surface(bottom_surface_oid)

    with project.edit(object_oid) as data_object:
        matched_primitives = spatial_filter.filter(data_object)

        if hasattr(data_object, "point_colours"):
            data_object.point_colours[matched_primitives] = [220, 220, 0, 255]
        elif hasattr(data_object, "block_colours"):
            data_object.block_colours[matched_primitives] = [220, 220, 0, 255]
        else:
            # This should be unreachable, as the call to `filter()` should
            # have failed in this case.
            raise RuntimeError("Unable to filter object!")

if __name__ == "__main__":
    with Project() as main_project:
        with main_project.undo():
            main(main_project)

The following animation demonstrates this script in action:

All spatial filters support being combined with the | operator.

Combining multiple filters with a logical OR

The any_filters() function combines multiple spatial filters into a single filter that matches points or blocks matched by any of the sub-filters—essentially performing a logical OR across the group. A common use case is combining multiple solids—such as ore bodies within a container—into a single filter in a simple, convenient way.

The following snippet demonstrates combining three spatial filters into a single filter using any_filters():

inside_any_filter = spatial_filters.any_filters(
    spatial_filters.inside_solid(surface_a),
    spatial_filters.inside_solid(surface_b),
    spatial_filters.inside_solid(surface_c),
)

The same result could be achieved by combining the filters using the| operator, as follows:

inside_any_filter = (
    spatial_filters.inside_solid(surface_a) |
    spatial_filters.inside_solid(surface_b) |
    spatial_filters.inside_solid(surface_c)
)

As discussed previously with all_filters(), the main advantage of any_filters() is that it allows the sub-filters to be created dynamically.

The animation below demonstrates a script in action using any_filters(). Blocks located inside any of the three yellow cubes are coloured yellow; all other blocks are coloured cyan:

Inverting filter classifications

The classification provided by a spatial filter can be inverted to select points or blocks that are not matched by the original filter. You can achieve this either by inverting the Boolean array returned by the filter or by creating a new inverted filter using the SpatialFilter.invert() method. Inverting the array is usually faster, while the invert() method is useful when combining filters programmatically.

Inverting the Boolean array

The filter() function of a spatial filter returns a NumPy Boolean array, with one element for each point or block in the input object: True indicates a match, and False indicates no match. Because this is a NumPy array A multidimensional, homogeneous array of values provided by the NumPy library, allowing efficient storage and computation on numerical data in Python., it can be quickly inverted using the tilde operator (~). The following script uses a filter to colour points above a surface yellow, then inverts the filter to colour points below the surface cyan:

from mapteksdk.data import Surface, spatial_filters
from mapteksdk.operations import object_pick
from mapteksdk.project import Project

def main(project: Project):
    surface_oid = object_pick(
        object_types=(Surface), label="Pick surface to filter blocks above"
    )
    object_oid = object_pick(label="Pick object to filter")

    spatial_filter = spatial_filters.above_surface(surface_oid)

    with project.edit(object_oid) as data_object:
        matched_primitives = spatial_filter.filter(data_object)

        if hasattr(data_object, "point_colours"):
            data_object.point_colours[matched_primitives] = [220, 220, 0, 255]
            data_object.point_colours[~matched_primitives] = [0, 220, 220, 255]
        elif hasattr(data_object, "block_colours"):
            data_object.block_colours[matched_primitives] = [220, 220, 0, 255]
            data_object.block_colours[~matched_primitives] = [0, 220, 220, 255]
        else:
            # This should be unreachable, as the call to `filter()` should
            # have failed in this case.
            raise RuntimeError("Unable to filter object!")

if __name__ == "__main__":
    with Project() as main_project:
        with main_project.undo():
            main(main_project)

The following animation shows this script in action:

Inverting a spatial filter

Alternatively, the SpatialFilter.invert() method can be used to create a new spatial filter that matches any point or block not matched by the original filter. In many cases, it is more efficient to invert the Boolean array returned by filter() using NumPy, because using invert() creates a new filter, which effectively runs the original filter twice and can double the script’s runtime.

The invert() method is most useful when a filter needs to be inverted before being combined with other filters using & or |. Note that inverting a filter is not the same as using the opposite filter. For example, the following filter:

spatial_filters.above_surface(surface_id).invert()

is not equivalent to this filter:

spatial_filters.below_surface(surface_id)

The correct equivalent is as follows:

spatial_filters.below_surface(surface_id, include_on=True) | spatial_filters.outside_surface_extent(surface_id)

This is because it inverts both the side of the surface and the inclusion of points on or outside the surface.