Module `gplately.ptt.utils.proximity_query`

Efficient distance queries involving many points tested against geometries.

Uses a point spatial tree to improve efficiency over single point/geometry queries. Most efficient when points are relatively uniformly spaced points.

EXAMPLE 1: Find the closest geometry(s) to each point in a sequence of points.

    import proximity_query
    import math

    # A list of 'pygplates.PointOnSphere' points.
    points = [...]

    # Some geometry features (eg, coastlines).
    geometry_feature_collection = pygplates.FeatureCollection('geometries.gpml')

    # Look for features within 90 degrees of each point.
    distance_threshold_radians = math.radians(90.0)

    # Extract the geometries from the features.
    geometries = []
    geometry_features = []
    for geometry_feature in geometry_feature_collection:
        geometries.append(geometry_feature.get_geometry())
        geometry_features.append(geometry_feature)


    # Find the closest geometry (feature) to each point (within threshold distance).

    geometry_features_closest_to_points = proximity_query.find_closest_geometries_to_points(
            points,
            geometries,
            geometry_features,
            distance_threshold_radians = distance_threshold_radians)

    # Print name of closest feature to each point (if any).
    for point_index, closest_geometry_feature in enumerate(geometry_features_closest_to_points):
        if closest_geometry_feature is not None:
            distance, geometry_feature = closest_geometry_feature
            print('Closest feature to', points[point_index].to_lat_lon(), 'is', geometry_feature.get_name(),
                    'with distance', distance * pygplates.Earth.mean_radius_in_kms, 'kms')
        else:
            print('No features close to', points[point_index].to_lat_lon())

    #
    # Find all geometries (features) near each point (within threshold distance).
    #
    geometry_features_closest_to_points = proximity_query.find_closest_geometries_to_points(
            points,
            geometries,
            geometry_features,
            distance_threshold_radians = distance_threshold_radians,
            all_geometries=True)

    # Print names of the closest features to each point (if any).
    for point_index, geometry_feature_list in enumerate(geometry_features_closest_to_points):
        if geometry_feature_list:
            print('Closest features to', points[point_index].to_lat_lon(), 'are...')
            for distance, geometry_feature in geometry_feature_list:
                print('    ', geometry_feature.get_name(), 'with distance',
                        distance * pygplates.Earth.mean_radius_in_kms, 'kms')
        else:
            print('No features close to', points[point_index].to_lat_lon())

EXAMPLE 2: Find the closest point(s) to each geometry in a sequence of geometries.

    import proximity_query
    import math

    # Some geometry features (eg, coastlines).
    geometry_feature_collection = pygplates.FeatureCollection('geometries.gpml')

    # Extract the geometries from the features.
    geometries = []
    geometry_features = []
    for geometry_feature in geometry_feature_collection:
        geometries.append(geometry_feature.get_geometry())
        geometry_features.append(geometry_feature)

    # Some uniformly spaced lat/lon points.
    points = []
    lon_lat_points = []
    for lat in range(-90, 91):
        for lon in range(-180, 181):
            points.append(pygplates.PointOnSphere(lat, lon))
            lon_lat_points.append((lon, lat))

    # Look for points within 5 degrees of each geometry.
    distance_threshold_radians = math.radians(5.0)

    #
    # Find the closest point (lon, lat) to each geometry (within threshold distance).
    #
    lon_lat_points_closest_to_geometries = proximity_query.find_closest_points_to_geometries(
            geometries,
            points,
            lon_lat_points,
            distance_threshold_radians = distance_threshold_radians)

    # Print longitude/latitude of closest point to each geometry (if any).
    for geometry_index, closest_lon_lat_point in enumerate(lon_lat_points_closest_to_geometries):
        if closest_lon_lat_point is not None:
            distance, lon_lat_point = closest_lon_lat_point
            print('Closest point to geometry', geometry_features[geometry_index].get_name(), 'is',
                     lon_lat_point, 'with distance', distance * pygplates.Earth.mean_radius_in_kms, 'kms')
        else:
            print('No points close to geometry', geometry_features[geometry_index].get_name())

    #
    # Find all points (lon, lat) near each geometry (within threshold distance).
    #
    lon_lat_points_closest_to_geometries = proximity_query.find_closest_points_to_geometries(
            geometries,
            points,
            lon_lat_points,
            distance_threshold_radians = distance_threshold_radians,
            all_points=True)

    # Print longitude/latitude of the closest points to each geometry (if any).
    for geometry_index, lon_lat_point_list in enumerate(lon_lat_points_closest_to_geometries):
        if lon_lat_point_list:
            print('Closest points to geometry', geometry_features[geometry_index].get_name(), 'are...')
            for distance, lon_lat_point in lon_lat_point_list:
                print('    ', lon_lat_point, 'with distance', distance * pygplates.Earth.mean_radius_in_kms, 'kms')
        else:
            print('No points close to geometry', geometry_features[geometry_index].get_name())

Functions

def find_closest_geometries_to_points(points, geometries, geometry_proxies=None, distance_threshold_radians=None, return_closest_position=False, return_closest_index=False, geometries_are_solid=False, all_geometries=False, subdivision_depth=4)

Efficient point-to-geometry distance queries when there are many relatively uniformly spaced points to be tested against geometries.

Parameters

points : a sequence of 'pygplates.PointOnSphere': a sequence of points
geometries : a sequence of 'pygplates.GeometryOnSphere': a sequence of geometries
geometry_proxies : sequence of objects associated with 'geometries', optional: If not specified then the proxies default to the geometries themselves. These can be any object (such as the 'pygplates.Feature' that the geometry came from).
distance_threshold_radians : number, optional: Optional distance threshold in radians - threshold should be in the range [0,PI] if specified.
return_closest_position : bool, default=False: Whether to also return the closest point on each geometry - default is False.
return_closest_index : bool, default=False: Whether to also return the index of the closest point (for multi-points) or the index of the closest segment (for polylines and polygons) - default is False.
geometries_are_solid : bool, default=False: Whether the interiors of the geometries are solid or not - only applies to polygon geometries - default is False.
all_geometries : bool, default=False: Whether to find all geometries near each point (within threshold distance) or just the closest. Defaults to False (only returns closest geometry to each point).
subdivision_depth : number: The depth of the lat/lon quad tree used to speed up point-to-geometry distance queries. The lat/lon width of a leaf quad tree node is (90 / (2^subdivision_depth)) degrees. Generally the denser the 'points' the larger the depth should be. Setting this value too high causes unnecessary time to be spent generating a deep quad tree. Setting this value too low reduces the culling efficiency of the quad tree. However a value of 4 seems to work quite well for a uniform lat/lon spacing of 'points' of 1 degree and below without the cost of generating a deep quad tree. So most of the time the subdivision depth can be left at its default value.

Returns

A list of geometry proxies associated with 'points'

The length of the returned list matches the length of 'points'. For each point in 'points', if the point is close to a geometry then that geometry's proxy (and its distance information) is stored (otherwise None is stored) at the same index (as the point) in the returned list. If 'all_geometries' is False then each item in returned list is a single geometry proxy (and its distance information) representing the closest geometry within threshold distance (or a single None). If 'all_geometries' is True then each item in returned list is a list of geometry proxies (and their distance informations) representing all geometries within threshold distance (or a single None). Above we mentioned "geometry proxy (and its distance information)". This is a tuple whose size depends on the values of 'return_closest_position' and 'return_closest_index' according to…

if return_closest_position and return_closest_index:
    geometry_proxy_to_point = (distance, closest_position, closest_index, geometry_proxy)
elif return_closest_position:
    geometry_proxy_to_point = (distance, closest_position, geometry_proxy)
elif return_closest_index:
    geometry_proxy_to_point = (distance, closest_index, geometry_proxy)
else:
    geometry_proxy_to_point = (distance, geometry_proxy)

The arguments 'distance_threshold_radians', 'return_closest_position', 'return_closest_index' and 'geometries_are_solid' are similar to those in pygplates.GeometryOnSphere.distance(). See http://www.gplates.org/docs/pygplates/generated/pygplates.GeometryOnSphere.html#pygplates.GeometryOnSphere.distance.

Raises ValueError if the lengths of 'geometries' and 'geometry_proxies' (if specified) do not match.

def find_closest_geometries_to_points_using_points_spatial_tree(points, spatial_tree_of_points, geometries, geometry_proxies=None, distance_threshold_radians=None, return_closest_position=False, return_closest_index=False, geometries_are_solid=False, all_geometries=False)

Same as 'find_closest_geometries_to_points()' except 'spatial_tree_of_points' is a 'points_spatial_tree.PointsSpatialTree' of 'points'.

This is useful when re-using a single 'points_spatial_tree.PointsSpatialTree'. For example, when using it both for point-in-polygon queries and minimum distance queries.

def find_closest_points_to_geometries(geometries, points, point_proxies=None, distance_threshold_radians=None, return_closest_position=False, return_closest_index=False, geometries_are_solid=False, all_points=False, subdivision_depth=4)

Efficient geometry-to-point distance queries when geometries are tested against many relatively uniformly spaced points.

Parameters

geometries : a sequence of 'pygplates.GeometryOnSphere': a sequence of geometries
points : a sequence of 'pygplates.PointOnSphere': a sequence of points
point_proxies : sequence of objects associated with 'points', optional: If not specified then the proxies default to the points themselves. These can be any object (such as a scalar value associated with the point).
distance_threshold_radians : number, optional: Optional distance threshold in radians - threshold should be in the range [0,PI] if specified.
return_closest_position : bool, default=False: Whether to also return the closest point on each geometry - default is False.
return_closest_index : bool, default=False: Whether to also return the index of the closest point (for multi-points) or the index of the closest segment (for polylines and polygons) - default is False.
geometries_are_solid : bool, default=False: Whether the interiors of the geometries are solid or not - only applies to polygon geometries - default is False.
all_points : bool, default=False: Whether to find all points near each geometry (within threshold distance) or just the closest. Defaults to False (only returns closest point to each geometry).
subdivision_depth : number: The depth of the lat/lon quad tree used to speed up point-to-geometry distance queries. The lat/lon width of a leaf quad tree node is (90 / (2^subdivision_depth)) degrees. Generally the denser the 'points' the larger the depth should be. Setting this value too high causes unnecessary time to be spent generating a deep quad tree. Setting this value too low reduces the culling efficiency of the quad tree. However a value of 4 seems to work quite well for a uniform lat/lon spacing of 'points' of 1 degree and below without the cost of generating a deep quad tree. So most of the time the subdivision depth can be left at its default value.

Returns

A list of point proxies associated with 'geometries'

The length of the returned list matches the length of 'geometries'. For each geometry in 'geometries', if the geometry is close to a point then that point's proxy (and its distance information) is stored (otherwise None is stored) at the same index (as the geometry) in the returned list. If 'all_points' is False then each item in returned list is a single point proxy (and its distance information) representing the closest point within threshold distance (or a single None). If 'all_points' is True then each item in returned list is a list of point proxies (and their distance informations) representing all points within threshold distance (or a single None). Above we mentioned "point proxy (and its distance information)". This is a tuple whose size depends on the values of 'return_closest_position' and 'return_closest_index' according to…

if return_closest_position and return_closest_index:
    point_proxy_to_point = (distance, closest_position, closest_index, point_proxy)
elif return_closest_position:
    point_proxy_to_point = (distance, closest_position, point_proxy)
elif return_closest_index:
    point_proxy_to_point = (distance, closest_index, point_proxy)
else:
    point_proxy_to_point = (distance, point_proxy)

Raises ValueError if the lengths of 'points' and 'point_proxies' (if specified) do not match.

def find_closest_points_to_geometries_using_points_spatial_tree(geometries, points, spatial_tree_of_points, point_proxies=None, distance_threshold_radians=None, return_closest_position=False, return_closest_index=False, geometries_are_solid=False, all_points=False)

Same as 'find_closest_points_to_geometries()' except 'spatial_tree_of_points' is a 'points_spatial_tree.PointsSpatialTree' of 'points'.

This is useful when re-using a single 'points_spatial_tree.PointsSpatialTree'. For example, when using it both for point-in-polygon queries and minimum distance queries.

def find_closest_points_to_geometry(geometry, points, point_proxies=None, distance_threshold_radians=None, return_closest_position=False, return_closest_index=False, geometry_is_solid=False, all_points=False, subdivision_depth=4)

Same as 'find_closest_points_to_geometries()' except with a single geometry instead of a list of geometries.

Returns a single "point proxy (and its distance information)" instead of a list (see 'find_closest_points_to_geometries()' for more information).

def find_closest_points_to_geometry_using_points_spatial_tree(geometry, points, spatial_tree_of_points, point_proxies=None, distance_threshold_radians=None, return_closest_position=False, return_closest_index=False, geometry_is_solid=False, all_points=False)

Same as 'find_closest_points_to_geometries_using_points_spatial_tree()' except with a single geometry instead of a list of geometries.

Returns a single "point proxy (and its distance information)" instead of a list.