Geo::Index - Geographic indexer
This document describes Geo::Index version 0.0.8
# Create and populate a geographic index use Geo::Index; @points = ( { lat => 1.0, lon => 2.0 }, { lat => -90.0, lon => 0.0, name => 'South Pole' }, { lat => 30.0, lon => -20.0, ele => 123.4 } ); $point = { lat=>10.0, lon=>20.0 }; $index = Geo::Index->new(); $index->IndexPoints( \@points ); $index->Index( $point ); $index->Index( [ 30, 40 ] ); # Search index %search_options = ( sort_results => 1, radius=>5_000_000 ); $results = $index->Search( [ -80, 20 ], \%search_options ); print "$$results[0]{name}\n"; # Prints 'South Pole' # Get all points in the southern hemisphere $results = $index->SearchByBounds( [ -180, -90, 180, 0 ] ); print "$$results[0]{name}\n"; # Also prints 'South Pole' ($closest) = $index->Closest( [ -80, 20 ] ); print "$$closest{name}\n"; # Also prints 'South Pole' ($closest) = $index->Closest( $points[1], { post_condition=>'NONE' } ); print "$$closest{name}\n"; # Also prints 'South Pole' ($farthest) = $index->Farthest( [ 90, 0 ] ); print "$$farthest{name}\n"; # Also prints 'South Pole' # Compute distance in meters between two points (using haversine formula) $m = $index->Distance( { lat=>51.507222, lon=>-0.1275 }, [ -6.2, 106.816667 ] ); printf("London to Jakarta: %i km\n", $m / 1000); $index->DistanceFrom( [ 90, 0 ] ); $m = $index->DistanceTo( $points[1] ); printf("Pole to pole: %i km\n", $m / 1000);
Geo::Index is a Perl module for creating in-memory geographic points indices. Once points have been indexed, fast searches can be run.
Efficient searches methods include Search(...) to get all points within a distance from a given point, SearchByBounds(...) to get all points in an given area, Closest(...) to get the closest points to a given point, and Farthest(...) to get the farthest points from a given point.
Search(...)
SearchByBounds(...)
Closest(...)
Farthest(...)
Additional methods are provided to compute distances between arbitrary points ( Distance(...), DistanceFrom(...), and DistanceTo(...) ) and to get the size in meters of one degree or the size in degrees of one meter at a given point (OneDegreeInMeters(...) and OneMeterInDegrees(...), respectively).
Distance(...)
DistanceFrom(...)
DistanceTo(...)
OneDegreeInMeters(...)
OneMeterInDegrees(...)
While by default computations are done for the Earth, other bodies can be used by supplying appropriates radii and circumferences to new(...).
new(...)
Geo::Index works with points on a spherical body. Points are hash references containing, at a minimum, lat and lon entries which give the point's position in degrees. Additional hash entries can be present and will be both ignored and preserved. The Index(...), IndexPoints(...), Search(...), Closest(...), Farthest(...), Distance(...), DistanceFrom(...), and DistanceTo(...) methods add additional entries in point hashes.
lat
lon
Index(...)
IndexPoints(...)
The hash entries used by Geo::Gpx are shown below. Apart from lat and lon these values are created by Geo::Gpx. Unless noted, these values may be read but should not be set, altered, or deleted.
lat - Point's latitude in degrees [ -90 .. 90 ]
lon - Point longitude in degrees [ -180 .. 180 )
These two values may be changed but the altered point should then be re-indexed using Index(...) before further searches are run.
data - The optional user data supplied when a point was created using the array shorthand. This contents of this field may be freely modified by the user. See Index(...) and IndexPoints(...), below.
data
lat_rad - The point's latitude in radians [ -pi/2 .. pi/2 ]
lat_rad
lon_rad - The point's longitude in radians [ -pi .. pi )
lon_rad
circumference - Circumference (in meters) of the circle of latitude that the point falls on. This is computed from the body's equatorial circumference assuming a spherical (not an oblate) body.
circumference
search_result_distance - Distance (in meters) of point from search point of previous search. The distance computation assumes a spherical body and is computed using a ruggedized version of the haversine formula. This value is only generated when Search(...) is called with the radius or sort_results option. See also Distance(...), DistanceFrom(...), and DistanceTo(...).
search_result_distance
radius
sort_results
antipode_distance - Distance (in meters) of point from search point's antipode as determined by a previous call to Farthest(...). This distance is computed using a ruggedized version of the haversine formula.
antipode_distance
As a convenience, most methods allow points to be specified using a shorthand notation [ lat, lon ] or [ lat, lon, data ]. Points given in this notation will be converted to hash-based points. If a point created using this notation is returned as a search result it will be as a reference to the hash constructed by Geo::Index and not as a reference to the original array. To access the data field of a point created using the shorthand notation use $$point{'data'} where $point is a search result point.
[ lat, lon ]
[ lat, lon, data ]
$$point{'data'}
$point
Any fields added to the indexed points by Geo::Index can be removed using Sweep(...) and Vacuum(...).
Sweep(...)
Vacuum(...)
$index = Geo::Index->new();
$index = Geo::Index->new()
Create a new empty index using default options: radius and circumferences are those of Earth, levels is set to 20 (~40 m index resolution).
levels
$index = Geo::Index->new( \@points );
Create a new index using default options and populate it with the given points.
The points in the array can be in either hash or array notation.
$index = Geo::Index->new( \%options );
Create a new empty index using specified options.
$index = Geo::Index->new( \@points, \%options );
Create a new index using specified options and populate it with the given points.
The options hash:
When a Geo::Index object is created, one can specify various options to fine-tune its behavior. The default values are suitable for a high-resolution index of Earth.
Average planetary radius (in meters). (default: 6371230)
If a radius is specified but polar_circumference or equatorial_circumference are not given then they will be calculated from the radius ( 2 * pi * radius )
polar_circumference
equatorial_circumference
Polar (meridional) circumference of the object the points lie on (in meters). (default: 40007863)
Circumference at the equator of the object the points lie on (in meters). (default: 40075017)
Depth of index. (valid: >0, <31; default: 20)
Note that the levels parameter specifies the number of non-full-globe index levels to generate and NOT the deepest index level. (Level -1, covering the entire globe, is always generated) For example, setting levels to 20 generates indices at levels 0 through 19 (plus level -1).
A summary of typical tile levels is shown below. To choose a value for the levels option using the table add 1 to the 'Level' shown for the desired maximum level of detail. The 'Grid' column shows the north-south size of each tile in meters at a each level. The 'Size' column shows the initial amount of RAM needed for an indexed set of 1 million random points using numeric keys on a 64-bit system when that level is the most detailed one (sizes may grow moderately once searches are run).
Level Grid Size Level Grid Size ----- --------- --------------- ----- ------- ------- -1 ~40000 km (entire planet) 12 ~5 km ~2.4 GB 0 ~20000 km 13 ~2.5 km ~2.7 GB 1 ~10000 km ~1.0 GB 14 ~1.2 km ~3.1 GB 2 ~5000 km ~1.0 GB 15 ~600 m ~3.3 GB 3 ~2500 km ~1.1 GB 16 ~300 m ~3.6 GB 4 ~1250 km ~1.2 GB 17 ~150 m ~3.8 GB 5 ~625 km ~1.3 GB 18 ~75 m ~4.1 GB 6 ~315 km ~1.4 GB 19 ~40 m ~4.4 GB 7 ~155 km ~1.5 GB 20 ~20 m ~4.6 GB 8 ~80 km ~1.6 GB 21 ~10 m ~4.9 GB 9 ~40 km ~1.7 GB 22 ~5 m ~5.1 GB 10 ~20 km ~1.9 GB 23 ~2 m ~5.4 GB 11 ~10 km ~2.1 GB 24 ~1 m ~5.6 GB
For reference, the memory usage of the array of 1 million random, unindexed points is about 440 MB, growing to about 540 MB with index use (about 100 bytes per point); the former amount is included in the index memory usage shown above.
function_type
Choose the type of low-level functions to use. (default: 'double' if available, 'perl' otherwise)
double
perl
Geo::Index will attempt to use compiled C code to speed up certain calculations. If the compilation fails (or was blocked by the user) then equivalent (but slower) Perl code will be used.
This option can be used to explicitly request the type of code to use. When set to 'float' then compiled C code using single-precision floating point will be requested. When set to 'double' then compiled C code using double-precision floating point will be requested. When set to 'perl' then Perl code will be used. If compiled code is unavailable then Perl code will be used regardless of what was requested.
float
Perl natively uses double-precision floating point. On modern hardware double-precision is slightly faster than single-precision. On certain platforms, however, it may be preferable to use single-precision instead of double-precision floating point. When needed, using single-precision should not be an issue since the minor errors introduced from loss of precision are drowned out by the errors inherent in the haversine function that is used for distance calculations.
$index->IndexPoints( \@points );
Add points in list to the index
If a point is added that already exists in the index and its position has changed then the existing index entry will be deleted and the point will be indexed again. If its position has not changed then no action will be taken.
@points
The points to add to the index
Each point in the list is either a reference to a hash containing at a minimum a lat and a lon value (both in degrees) or a reference to an array giving the point. See the Points section above for details.
$index->Index( \%point );
$index->Index( \@point );
Add a single point to the index
If the point being added already exists in the index and its position has changed then the existing index entry will be deleted and the point will be indexed again. If its position has not changed then no action will be taken.
%point or @point
%point
@point
The point to add to the index
This can be either a hash containing at a minimum a lat and a lon value (both in degrees) or an array giving the point. See the Points section above for details.
$index->Unindex( \%point );
Remove specified point from index
This method will remove the point from the index but will not destroy the actual point.
The point to remove from the index.
Note that this must be a reference to the actual point to remove and not to a copy of it. Simply specifying a point's location as a new point hash will not work.
Added in 0.0.4 to replace the functionally identical (and now deprecated) DeletePointIndex(...).
DeletePointIndex(...)
@results = $index->Search( \%point, \%options );
$results = $index->Search( \%point, \%options );
Search index for points near a given point
The point to search near
This is either a reference to a hash containing at a minimum a lat and a lon value (both in degrees) or a reference to an array giving the point. See the Points section above for details.
%options
The parameters for the search (all are optional):
Note that except for radius, none of the below options have any effect when quick_results is specified.
quick_results
Only return results within this distance (in meters) from search point.
If no radius is specified or the radius is set to Geo::Index::ALL then all points in the index will be returned.
Geo::Index::ALL
When quick_results is specified then all points within the specified radius will be returned (additional points outside the radius may also be returned).
Sort results by distance from search point before filtering and returning them.
max_results
Return at most this many results.
Unless sorting is also requested these are not guaranteed to be the closest results to the search point.
pre_condition
Reference to additional user-supplied code to determine whether each point should be included in the results.
This code is run before the distance from the search point to the result point has been calculated.
See the Condition functions section below for syntax.
post_condition
This code is run after the distance from the search point to the result point has been calculated.
user_data
Arbitrary user-supplied data that is passed to the condition functions.
This can be used to allow the function access to additional data structures.
Return the raw preliminary results. (no result filtering is done)
Normally when a search is performed the raw results have their distances computed and are filtered against the search radius and the condition functions. This process can be very slow when there are large numbers of points in a result set. Specifying this option skips those steps and can make searches several orders of magnitude faster.
When this option is active then instead of returning a list of points the Search(...) method will return a list of lists of points (some of which may be undef). If a radius was specified then the results returned will contain all the points within that radius. Additional points outside of the search radius will likely be returned as well.
undef
When iterating over the arrays be sure to check whether a list element is undef before trying to deference it.
An example of returned quick results (in scalar context); POINTs are references to different points:
[ [ POINT, POINT, POINT ], [ POINT, POINT ], undef, [ POINT, POINT ] ]
To be clear, when this option is active rough radius limiting is done but there is no filtering done, no distances are computed, and no sorting is performed.
See the Performance section below for a discussion of this option and when to use it.
tile_adjust
Manual adjustment for tile size (signed integer, default is 0)
0
Values < 0 use smaller tiles, values > 0 use larger tiles. Each increment of 1 doubles or halves the tile size used. For example, set to -1 to use tiles half the normal size in each direction.
< 0
> 0
1
-1
This option can be a bit counter-intuitive. Although using smaller tiles will result in fewer points that need to be checked or returned it will also lead to a larger number of tiles that need to be processed. This can slow things down under some circumstances. Similarly using larger tiles results in more points spread over fewer tiles. What adjustment (if any) will result in the highest performance is highly dependant on both the search radius and on the number and distribution of the indexed points. If you adjust this value be sure to benchmark your application using a real dataset and the parameters (both typical and worst-case) that you expect to use.
Return value
In list context the return value is a list of references to the points found or an empty array if none were found.
In scalar context the return value is a reference to the aforementioned list or undef if no results were found.
If either the sort_results or radius options were specified in the search options then for each point in the results the distance in meters from it to the search point will be stored in the search_result_distance entry in the result point's hash. It can be retrieved using e.g. $meters = $$point{search_result_distance};
$meters = $$point{search_result_distance};
See above section for the results returned when the quick_results option is active.
@results = $index->SearchByBounds( \@bounds, \%options );
@results = $index->SearchByBounds( \%bounds, \%options );
$results = $index->SearchByBounds( \@bounds, \%options );
$results = $index->SearchByBounds( \%bounds, \%options );
Search index for points within a given bounding box
The points returned are those that lie between the specified latitudes and longitudes. The four corners form a rectangle only when using certain map projections such as equirectangular or mercator (including pseudo-mercator a.k.a. web mercator as used by slippy maps). If you are using a projection that does not have horizontal lines of latitude and vertical lines of longitude and you want the results to lie within and/or fill a rectangle on your map then your code will need to perform the appropriate bounds adjustment and point filtering itself.
@bounds or %bounds
@bounds
%bounds
The search boundaries
This can be specified either as a list or a hash. Any of the following are acceptable:
( w_val, s_val, e_val, n_val )
( 'south' => s_val, 'north' => n_val, 'west' => w_val, 'east' => e_val )
( 's' => s_val, 'n' => n_val, 'w' => w_val, 'e' => e_val )
( 'S' => s_val, 'N' => n_val, 'W' => w_val, 'E' => e_val )
s_val and n_val are the south and north latitudes of the bounding box and w_val and e_val are its west and east longitudes (all values are in degrees). For the list form the order matches that used by PostGIS, shapefiles, and GeoJSON but be aware that the order of the fields is not standardized.
s_val
n_val
w_val
e_val
Note that none of the below options have any effect when quick_results is specified.
condition
Note that unlike with the other search methods there is only a single condition function. Instead of the $_search_point, the second parameter to the condition function is a reference to the bounding box in list form (as described above for @bounds). Lastly, since "distance from search point" makes no sense in the context of a bounding box, none is provided to the function.
$_search_point
Arbitrary user-supplied data that is passed to the condition function.
Normally when a search is performed the raw results are filtered to ensure that they lie within the bounds and the condition function is applied to each result point. This process can be slow when there are large numbers of points in a result set. Specifying this option skips those steps and can make searches faster.
When this option is active then instead of returning a list of points the SearchByBounds(...) method will return a list of lists of points (some of which may be undef). The results are guaranteed to contain all points within the requested bounds but additional points outside of the requested bounds will likely be returned as well.
To be clear, when this option is active rough bounds limiting is done but there is no filtering done and no bound checks are actually performed.
@results = $index->Closest( \%point, $number_of_points_desired, \%options );
$results = $index->Closest( \%point, $number_of_points_desired, \%options );
Find the point or points closest to a given point
Note that if you want to find the closest points within a given radius it may be faster to use Search(...) instead. See the Performance section below for more details.
$number_of_points_desired
The number of points that should be returned.
Set to 0 to not restrict the number of points returned or set it > 0 to set the maximum number of points to return.
If omitted then this will default to 1.
If no radius is specified or the radius is set to Geo::Index::ALL then all points in the index may potentially be returned.
Sort results by distance from point
By default points returned are sorted by distance. Set this to 0 to not sort the returned points.
Although sorting is not mandatory, performing it is strongly recommended since otherwise the set of points returned are not guaranteed to be the closest.
By default, a post_condition function that filters out the search point is used. To remove this default function either specify a new one or set post_condition to "NONE".
NONE
For each point in the results the distance in meters from it to the search point will be stored in the search_result_distance entry in the result point's hash. It can be retrieved using e.g. $meters = $$point{search_result_distance};
@results = $index->Farthest( \%point, $number_of_points_desired, \%options );
$results = $index->Farthest( \%point, $number_of_points_desired, \%options );
Find the point or points farthest from a given point
In other words, find the points closest to a given point's antipode.
The point to search relative to
Set to 0 to not restrict the number of points returned or set it >0 to set the maximum number of points to return.
Although sorting is not mandatory, performing it is strongly recommended since otherwise the set of points returned are not guaranteed to be the farthest.
By default, a post_condition function that filters out the search point is used. To remove this default function either specify a new one, set a value for user_data, or set post_condition to "NONE".
If the default post_condition is active and no user_data value has been provided by the caller then this is set to the actual (non-antipodal) search point.
For each point in the results the distance in meters from it to the search point will be stored in the search_result_distance entry in the result point's hash. In addition, the distance from a result point to the search point's antipode will be stored in the antipode_distance entry. These can be retrieved using e.g.:
$meters_from_search_point = $$point{search_result_distance}; $meters_to_antipodal_point = $$point{antipode_distance};
$meters = $index->Distance( \%point_1, \%point_2 );
$meters = $index->Distance( \@point_1, \@point_2 );
$meters = $index->Distance( \%point_1, \@point_2 );
$meters = $index->Distance( \@point_1, \%point_2 );
Returns the distance in meters between two points
The haversine function is used to compute the distance. As this assumes a spherical body the distances returned may show errors. Using the default options, these errors are up to 0.056% (north - south) or 1.12% (east - west). Such errors typically start becoming significant at distances over 20 km.
%point_1 or @point_1, %point_2 or @point_2
%point_1
@point_1
%point_2
@point_2
The points to measure the distance between
These can be either hashes containing at a minimum a lat and a lon value (both in degrees) or arrays giving each point. See the Points section above for details.
$meters = $index->DistanceFrom( \%point_1 );
$meters = $index->DistanceFrom( \@point_1 );
Set an initial point to measure distances from
Note that any call to Distance(...) and some calls to Search(...) (those using the radius or sort_results options) will overwrite the initial point set with this method.
%point_1 or @point_1
The point to measure distances from
$meters = $index->DistanceTo( \%point_2 );
$meters = $index->DistanceTo( \@point_2 );
Returns the distance in meters between the specified point and the one set earlier with DistanceFrom(...).
%point_2 or @point_2
The point to measure distances to
%configuration = $index->GetConfiguration( );
Returns the running configuration of the Geo::Index object.
See also GetStatistics(...) and examples/show_configuration.pl
GetStatistics(...)
examples/show_configuration.pl
The return value is a hash with the following entries:
key_type - The key type in use:
key_type
'text' for text keys (e.g. '12:345,6789')
text
12:345,6789
'numeric' for 64-bit numeric keys
numeric
'packed' for 64-bit numeric keys packed into an 8-byte string
packed
supported_key_types - The types of keys that can be used
supported_key_types
Value is a reference to a list of supported key types (as given above).
code_type - The type of low-level code in use:
code_type
'perl' for Perl functions
'float' for C functions mostly using float values.
'double' for C functions mostly using double values.
supported_code_types - The types of low-level code that can be used
supported_code_types
Value is a reference to a list of supported code types (as given above).
levels - Number of levels in index (excluding the global level)
planetary_radius - Average planetary radius (in meters)
planetary_radius
polar_circumference - Polar circumference (in meters)
equatorial_circumference - Equatorial circumference (in meters)
size - Number of points currently indexed
size
tile_meters - Size in meters (at the equator) of each tile the at most-detailed level of index
tile_meters
@stats = $index->GetStatistics( );
Returns statistics regarding the Geo::Index object.
See also GetConfiguration(...) and examples/show_configuration.pl
GetConfiguration(...)
The return value is a list with one entry per level. Each list entry is a hash reference giving statistics for a single level of the index and contains the following entries:
level - The level number the statistics are for
level
points - Total number of points indexed in this level
points
tiles - Number of tiles containing at least one point
tiles
min_tile_points - Minimum number of points in a non-empty tile
min_tile_points
max_tile_points - Maximum number of points in a non-empty tile
max_tile_points
avg_tile_points - Average number of points in a non-empty tile
avg_tile_points
$index->Sweep( );
$index->Sweep( \%point );
$index->Sweep( \@points );
$index->Sweep( undef, \@extra_keys );
$index->Sweep( \%point, \@extra_keys );
$index->Sweep( \@points, \@extra_keys );
Remove data generated by searches from some or all points
The fields that will be removed are search_result_distance and antipode_distance.
Called on its own (with no point or points specified) this method will remove data generated by searches from all points.
A list of additional keys to remove can optionally be supplied. To request vacuuming of all points with additional keys specified, use undef instead of \%point or \@points.
\%point
\@points
See also Vacuum(...).
%point or @points
The point or a list of points to remove metadata from.
@extra_keys
List of additional keys to remove
$index->Vacuum( );
$index->Vacuum( \%point );
$index->Vacuum( \@points );
$index->Vacuum( undef, \@extra_keys );
$index->Vacuum( \%point, \@extra_keys );
$index->Vacuum( \@points, \@extra_keys );
Remove all data generated by Geo::Index from some or all points
The fields that will be removed are: lat_rad, lon_rad, circumference, search_result_distance, antipode_distance.
Called on its own (with no point or points specified) this method will remove all generated data from all points.
See also Sweep(...).
$count = $index->PointCount( );
Returns the number of points currently in index
@all_points = $index->AllPoints( );
$all_points = $index->AllPoints( );
Returns all points currently in index
In list context the return value is a list of references to all points in the index or an empty array if there are no points in the index.
In scalar context the return value is a reference to the aforementioned list or a reference to an empty list if there are no points in the index.
$index->OneMeterInDegrees( $latitude );
Returns length in degrees of one meter (N/S and E/W) at given latitude
Values are approximate. The diameters are those for an oblate spheroid but the math assumes a sphere. As one approaches the poles these values get heavily distorted; code that uses them needs to take this into account.
See also OneDegreeInMeters(...).
$latitude
The latitude in radians
An two-element list containing the width and height in meters of one degree at the given latitude:
( degrees_of_latitude, degrees_of_longitude )
Returns length in meters of one degree (N/S and E/W) at given latitude
See also OneMeterInDegrees(...).
An two-element list containing the width and height in degrees of one meter at the given latitude:
( north_south_meters, east_west_meters )
Geo::Index uses CamelCase for its method names.
For those who prefer using snake case, alternate method names are provided:
Method Alternate name ----------------- -------------------- Index index IndexPoints index_points Unindex unindex BuildPoints build_points AddValue add_value GetValue get_value Search search SearchByBounds search_by_bounds Closest closest Farthest farthest AllPoints all_points Distance distance DistanceFrom distance_from DistanceTo distance_to OneMeterInDegrees one_meter_in_degrees OneDegreeInMeters one_degree_in_meters GetConfiguration get_configuration GetStatistics get_statistics Sweep sweep Vacuum vacuum PointCount point_count
The Search(...), SearchByBounds(...), Closest(...), and Farthest(...) methods allow a user-supplied condition function to filter potential results.
If present, these condition functions are called for each potential search result. They should be idempotent* and could potentially be called multiple times for a given point. The code should return TRUE (e.g. 1) if a potential point should be included in the results or FALSE (e.g. 0 or undef) if the point should be excluded.
For Search(...), Closest(...), and Farthest(...), the pre_condition function runs before the distance to the result point has been calculated and the post_condition function runs after it has been calculated. For SearchByBounds(...) no distances are calculated and the function is simply called once per point.
* Functions can set outside values provided they do not affect any values used internally by Search(...) and so long as those outside values have no effect on the condition's outcome. Such behavior is, of course, frowned upon.
The parameters to the condition function are, in order:
$_result_point
Reference to the potential search result being checked
Reference to the point at the center of the search
For SearchByBounds(...) this is instead the bounding box: [ west, south, east, north ]
[ west, south, east, north ]
$user_data
Arbitrary user-supplied data
For example, the options set in the following code allows all points in the results except for the one named 'Point Nada':
$options{pre_condition} = sub { my ( $_result_point, $_search_point, $user_data ) = @_; if ( $$_result_point{name} eq $user_data ) { return 0; # Exclude result } return 1; # Point is a valid search result }; $options{user_data} = "Point Nada";
To exclude the search point from the search results use:
$options{post_condition} = sub { my ( $_result_point, $_search_point, $user_data ) = @_; return ( $_result_point != $_search_point ); };
or more concisely
$options{post_condition} = sub { return $_[0] != $_[1]; };
In general, post_condition functions should be preferred since the overhead of the Perl function call is typically larger than that of the distance calculation. By checking the distance first, running the post_condition function might not be necessary.
Geo::Index is intended for stand-alone applications that need a way to quickly perform proximity searches on relatively small datasets (at most a few million points). Typical search speeds are three to five orders of magnitude faster than a linear search. For larger datasets and for applications running in a server environment using something like PostGIS is more appropriate.
Indexing speed is about 50,000 points per second when levels is 20. Search speeds are highly dependent on the data indexed and on search parameters but are typically in the neighborhood of a few thousand searches per second.
Memory usage tends to be rather high; for 1,000,000 points the index is ~3.2 GB for tightly clustered points or ~4.6 GB when spread evenly world-wide. The size of an index grows linearly with each added point at a rate of about 4 kB per point. When a point is first encountered whilst searching its size will increase by about 100 bytes (this only happens once per point).
Since performance is so dependant on data and usage, the the user is encouraged to test all available options while developing their application before choosing the one that works fastest. The examples/benchmark.pl script included with this module may be helpful for measuring this module's performance.
examples/benchmark.pl
Here are some guidelines for best results:
Requesting results as a list reference is faster than asking for a plain list.
That is, e.g., $results = Search(...); is faster than @results = Search(...);
$results = Search(...);
@results = Search(...);
Post-conditions are faster than pre-conditions.
Benchmarking has shown that the cost of the Perl function call is higher than that of the distance-related code. Thus there is probably no reason to use pre-conditions. Put concisely,
$results = $index->Search( $point, { post_condition => $code_ref } );
is faster than
$results = $index->Search( $point, { pre_condition => $code_ref } );
Choose an appropriate value for levels when creating the index
The Search(...) method has best performance when the size of the most detailed level of the index has a smaller physical size than the radius of a typical search. For example, if your searches are typically for points within 100 meters then an index with levels should be set to at least 18 (~75 meters at the equator) to yield best results; if typical searches have 10 meter radius then levels should be 22.
The Closest(...) method works best when the most detailed level of the index contains a single point per tile and search points lie close to potential result points.
To help tune the levels value, the GetConfiguration( ) method can be used to find out the physical size of the most detailed level along with statistics on the number of points per index tile.
GetConfiguration( )
Use the quick_results option when possible.
Filtering points and combining them into a single, flat list can be very expensive. Many applications can tolerate getting additional points beyond those matching the search criteria. An example of this is drawing points on a map; if points are clipped to the visible area when they are drawn it may not matter if some of them lie outside of it.
Use Search(...) instead of Closest(...) when you have a search radius.
The Closest(...) function is most efficient when no search radius is specified or when result points lie very close to the search point. Closeness is relative to the tile size of the most detailed index level; for the default index depth (20), "very close" is roughly within about 100 meters.
20
When clipping results to a maximal radius it is typically much faster to use Search(...) with the sort_results and max_results options*.
For example, to find the closest $n points within distance $d of a point $p it is usually much faster to use
$n
$d
$p
%options = ( max_results => $n, radius => $d, sort_results => 1, post_condition => sub { return $_[0] != $_[1]; } ); $results = $index->Search( $p, \%options );
instead of
$results = $index->Closest( $p, $n { radius => $d } );
* The post_condition shown in the example omits the search point from the results and is needed to fully emulate the behavior of Closest(...).
Both Search(...) and SearchByBounds(...) are very fast since they can find the relevant index tiles in linear time. Since the time needed to filter the results is directly proportional to the number of points retrieved from the index, best performance occurs when the size of the most detailed tiles is smaller than that of the typical search radius or search bounds.
Searches run using Closest(...) are done starting from the most detailed level and work upwards. Best performance occurs when a result is found in the first few iterations. If the first iteration that finds points yields a large number of points then performance will suffer since the distance to each of these points will need to be measured to find the closest. For similar reasons, requesting a large number of closest points in a single call will also impact performance. The Farthest(...) method is largely a wrapper for Closest(...) and thus exhibits similar behavior.
Some functions within Geo::Index have optional implementations written in C. If these are active (by default they are whenever possible) searches typically run 25% to 50% faster.
Whenever possible Geo::Index uses numeric index keys. Compared to text index keys, numeric keys improve performance with about 30% faster speed and about 50% smaller index memory usage. The downside to numeric keys is that they are less legible to humans while debugging. (Whether numeric or text keys are used can be changed by setting the appropriate value at the top of Geo/Index.pm)
Geo/Index.pm
Typical benchmark results run on a modern workstation using numeric keys and double-precision C low-level code with the index containing 1,000,000 points are as follows:
Points can be added to an index at the rate of about 50,000 per second.
Typical searches returning values run at about 25,000 to 50,000 searches per second. Worst-case performance is under 50 searches per second and searches returning no results run at over 100,000 searches per second. The overhead of traversing the results is fairly negligable.
Quick searches run at 120,000 to 150,000 searches per second. Actually doing anything with the results slows things down a lot. Including traversal of the results, a typical quick search runs at 40,000 to 100,000 searches per second with the worst-case being about 80 searches per second.
If distances to the result points are not needed, quick searches are typically about 75% faster than normal ones albeit with about 5 times as many results being returned.
For the SearchByBounds(...) method run time correlates with the size of the bounding box with smaller bounding boxes typically yielding faster run times.
A fairly typical search yielding about 50 results runs at about 10,000 searches per second in normal mode and about 30,000 searches per second in quick mode. A nearly worst case example is a search returning 100,000 points; this will run at about 5 searches per second in normal mode or about 8,000 searches per second in quick mode.
For the Closest(...) method the highest performance is seen when there are result points close to the search point. Search speeds for the single closest point are typically in excess of 20,000 per second for close-by results or about 8,000 per second when results are far away. Worst case speeds of about 1,000 searches per second occur when all indexed points are in the hemisphere opposite the search point.
For the Farthest(...) method the highest performance is seen when there are result points nearly antipodal to the search point. Search speeds for the single farthest point are typically in excess of 20,000 per second when nearly-antipodal points exist. Worst case speeds of about 1,000 searches per second occur when all indexed points are in the same hemisphere as the search point.
Note that the numbers above are approximate and are highly dependant on the data being searched, the type of search being run, and on the number of results returned. Actual searches may be an order of magnitude or more slower.
A sample benchmark run can be found in examples/benchmark.txt To run the benchmarks yourself you can run examples/benchmark.pl It needs the Devel::Size and Time::HiRes modules installed and a single run takes about 8 minutes.
examples/benchmark.txt
Since Perl constants cannot be changed from the commandline you will need to edit the Geo/Index.pm to force the use of numeric keys, packed numeric keys, or text keys. This can be done by uncommenting the appropriate lines at the head of the file (look for USE_NUMERIC_KEYS and USE_PACKED_KEYS). When running benchmark.pl, the various other options can be found at the top of the script. When writing your own programs you can switch between the Perl and C single- or double-precision low-level code by using the function_type option when calling new(...).
USE_NUMERIC_KEYS
USE_PACKED_KEYS
benchmark.pl
The high cost of constructing and traversing the results seems inherent to the Perl language and there does not seem to be any way to avoid it. The is some potential for optimization though:
The pre_condition and post_condition function calls might be sped up by assigning them to function handles (much as is done with HaversineDistance, etc.) instead of making the calls by dereferencing the variables.
HaversineDistance
Performance could potentially be improved by splitting the current combined index into individual indices for each level. Having smaller keys and indices should result in higher performance but the additional layer of indirection could slow things down in some circumstances.
Improvements might be possible to the performance of Closest( n, ...) where n>1 and per-point distances are not needed by the caller. At each iteration of the algorithm the previously-used search radius gives the maximal distance to all points already found, obviating the need to calculate every point's distance. Only points in the final level of the search would need to have their distances calculated. The downside to this method is that the point distances would not be available for all points in a result set (only for those found in the final search level).
Closest( n, ...)
n
A number of alternative datastructures were explored for the point index but benchmarking showed plain Perl hashes to be by far the most efficient. It is possible, though in my opinion unlikely, that a faster data structure choice exists that is suitable for use in this module.
A given index comprises sets of tiles at various zoom levels with each tile containing a list of the points that lie within it. The lowest level of the index covers the entire globe. Each higher index level contains twice as many tiles in each direction. At each zoom level points are linearly mapped to grid tiles based on their latitudes and longitudes using an equirectangular projection. This is fairly analogous to how typical web slippy maps are organized (though they use a pseudo-mercator projection).
As one approaches the poles the tiles become increasingly distorted with the area (in square meters) covered by each tile becoming progressively smaller. The distance in meters for one degree of longitude gets smaller as one moves away from the equator. The distance in meters for one degree of latitude, however, remains constant at all latitudes.
Each tile has a name that gives its zoom level and position. These names are used as keys into a Perl hash allowing the quick retrieval of the points that lie within a given tile. The various search methods are designed to efficiently pull points from this index and filter them in various ways. The format used for the keys is described in the Tile naming section below.
Additional datastructures (e.g. the list of all points in the index) are also present but knowing their details is not needed to understand how the index functions. In the descriptions below, some minor (but often critical) details have been omitted and some simplifications have been made; these details (mainly edge cases) are discussed in the code comments.
When a point is added to the index it is stored multiple times in the index hash, once for each level. This is done as follows:
The point's latitude and longitude are converted to integers. This is done using a simple linear mapping. In pseudo-code, the equations used are:
max_size = 2**levels integer_latitude = floor( ( latitude + 90.0 ) * max_size / 180.0 ) integer_latitude = max_size - 1 if (integer_latitude == max_size) integer_longitude = floor( ( longitude + 180.0 ) * max_size / 360.0 ) % max_size
The values for latitude and longitude are in degrees and levels is the number of levels in the index (not counting the global one).
latitude
longitude
Each index level is looped through from the index's maximum level to 0. At each level, the key (comprised of level, integer_latitide, and integer_longitide, see also below) is used to retrieve the corresponding value from the index hash. This value is a reference to the list of points that lie within the grid tile named by the key. The point being indexed is added to the retrieved list. If there is no list stored in the index for the current key then a new list is created and added. As a special case, all points adjacent to the poles (that is points with integer latitudes of 0 or max_size - 1) use the longitude ALL in their keys.
integer_latitide
integer_longitide
max_size - 1
ALL
Once the point has been added, the integer latitudes and longitudes as well as the max_size are shifted right by one bit in preparation for the the next level.
max_size
Once a the point has been added to the index at each level, the point is added to the global index entry using the key ALL, ALL, ALL. (All indexed points can be found under this key.)
The Search(...) method is typically used to find all points lying within a given radius of a search point. Two steps are performed by this method: retrieval of preliminary results and filtering of the results based on the search criteria.
If no search radius was specified, if a global search was requested, or if the search radius covers more than half the globe then the preliminary results are all points in the index. Otherwise, the preliminary results are gathered as follows:
The appropriate tile zoom level to use is determined using:
shift = ceil( log2( search_radius / half_circumference ) ) level = max_level - shift
This results in the smallest level that could potentially contain all result points within a single tile.
The search radius (in meters) is converted to two search angular radii, one for latitude and one for longitude. This is done since the number of meters per degree longitude decreases as one approaches the poles. Thus the north-south (latitude) search radius remains constant at all latitudes but the east-west (longitude) search radius increases as one nears the poles.
Each extreme is converted to an integer and shifted right by the determined shift, The preliminary results are retrieved from the index by iterating over the keys for the computed level, bounded by the integer extrema. This typically, but not always, results in a list of pointers to four tiles' points.
shift
If the quick_results option is active then this preliminary list of lists of points is returned. If not then the points are filtered to only include those matching the search criteria. The filtered points are optionally sorted and then returned. Note that when large numbers of points have been found this filtering can be very slow; see Performance above for details.
The Closest(...) and Farthest(...) methods find the points closest to (or farthest from) a search point. The Closest(...) method works as follows:
The search starts at the most detailed level of the index and proceeds to the least detailed (0). At each level, the grid tile that the search point lies in along with the three closest grid squares are identified. The method used for selecting the adjacent tiles is to look at the least-significant bits of the integer position at the previous (more detailed) level. A 1 bit for the latitude selects tiles to the north, a 0 bit the ones to the south. Likewise a 1 for the longitude selects the ones east and a 0 the ones west.
Now that the four tiles have been identified, the largest radius from the search point to the tile edges is determined. The distance from the search point to each point within the four tiles is measured. If the point is within the radius computed and it passes any pre- or post-condition tests it is added to the results list. To speed up processing, points that have already been rejected along with the distances so far measured are cached. As a convenience, by default a filter is applied that omits the search point from the results.
Once all points within a level's four chosen tiles have been gathered a check is done to see whether at least the requested number of points have been found. If they have then the loop ends, if not then the next (less-detailed) level is processed.
By default, the results are then sorted by distance which ensures that the closest results are earliest in the list. This is necessary since although the nature of the algorithm tends to place closer points earlier in the results there is no inherent order to the points added from a particular index level. Lastly, the requested number of result points is returned.
The Farthest(...) method is largely implemented as a wrapper for Closest(...). It functions by finding the closest points to the search point's antipode.
The SearchByBounds(...) method works much the same as Search(...) method. Instead of computing extrema of a search circle, those of the supplied bounding box are used. The tile level used is max( latitude_level, longitude_level ) where latitude_level and longitude_level are the most detailed levels that could potentially (given their extrema's angular distances) contain their respective extrema within a single tile in each direction. The remainder of the method is identical to that of Search(...) albeit with all distance-related code removed.
max( latitude_level, longitude_level )
latitude_level
longitude_level
As mentioned earlier, keys consist of a zoom level, a latitude, and a longitude. Each key uniquely names a given tile.
Zoom levels are either integers between 0 and the maximum zoom level or the special zoom level ALL (with the value -1) that covers the entire globe. Latitudes and longitudes are integers between 0 and one less than the maximum grid size for the level. The tiles immediately adjacent to the poles are treated differently. In these areas the coverage of each tile is quite small and the algorithm around the poles would normally be complex. To accommodate these issues, the special value ALL (with the value -1) is used for the longitude of the polar tiles (those areas with the lowest or highest latitude value for the key's level). All points lying in a polar region are assigned to that region's overarching ALL tile. At the global level all three components are set to ALL.
If Perl has been compiled with 64-bit support then each key is packed into a 64 bit integer. The level is stored in the upper 6 bits (bits 58 .. 63), the integer latitude in the next 29 bits (bits 29 .. 57), and the integer longitude in the low 29 bits (bits 0 .. 28). To represent the ALL value all bits in the relevant field are set to 1. Note that even on 32-bit systems Perl is often compiled with 64-bit support.
If Perl does not have 64-bit support then a different format is used. In most places within Geo::Index, keys are stored as three-element array references. The first field contains the level, the second the integer latitude and the third the integer longitude. If present, ALL values are stored as-is as their integer value (-1). For accessing the index, keys are converted to strings with the format "level:latitude,longitude" with the literal string "ALL" being used for ALL values.
:
,
Each index object is a hash containing a number of entries These are:
$self->{index} - The points index
$self->{index}
Entry is a hash reference. Keys are tile names (as discussed above), values are lists of point references.
$self->{indices} - Indices used for each point
$self->{indices}
Entry is a hash reference. Keys are point references, values are lists of tile names.
$self->{positions} - Each point's position when indexed
$self->{positions}
Entry is a hash reference. Keys are point references, values are two-element lists giving each point's latitude and longitude at the time it was indexed.
$self->{levels} - Number of levels in the index (excluding the global level)
$self->{levels}
$self->{max_level} - The highest-resolution level number (i.e. levels - 1)
$self->{max_level}
$self->{max_size} - Number of grid tiles in each direction at most detailed level of index
$self->{max_size}
$self->{planetary_radius} - The planetary radius used by the index (in meters)
$self->{planetary_radius}
$self->{polar_circumference} - The polar circumference used by the index (in meters)
$self->{polar_circumference}
$self->{equatorial_circumference} - The equatorial circumference used by the index (in meters)
$self->{equatorial_circumference}
This module is not believed to be thread-safe. In specific:
The SetUpDistance(...) function stores the first point's position in global variables.
SetUpDistance(...)
To fix this, DistanceFrom*(...) and DistanceTo*(...) would need to be removed plus SetUpDistance*(...) and HaversineDistance*(...) would need to be combined into a single 4-argument HaversineDistance*(...) function. Calling code would need to be modified as appropriate. In terms of performance, the overall cost of doing this is likely quite low.
DistanceFrom*(...)
DistanceTo*(...)
SetUpDistance*(...)
HaversineDistance*(...)
The search code stores distances computed for a specific search into the point datastructures. If multiple concurrent searches are run against a single index then distances computed by one search may overwrite those from another search. This can lead to inconsistent results.
To fix this a per-search distance hash would need to be maintained. This could have serious performance implications and would preclude returning the point distances within the point hashes. The distances could, however, be returned in an additional datastructure.
Adding and deleting points to/from the index is not atomic. Running e.g. a search while points are being added or deleted can lead to unpredictable behavior (up to and including the program crashing).
One could fix this by adding object-level locks:
Block concurrent calls to the Index(...) and Unindex(...)methods
Unindex(...)
Block calls to the Index(...) and Unindex(...) methods while searches are running
Block calls to Search(...) et al. when the Index(...) or Unindex(...) methods are active
Including the same point in multiple indices or searches at the same time could lead to interesting results.
As mentioned above, this is due to the storage of search result distances within the points and not within the index object. Each search that involves a given point will likely overwrite its search_result_distance value.
This could be encountered in a number of ways. For example, a search using a condition function that itself runs a search against the second index could be problematic. This could be encountered even when using a single index. For example, if code relies on the distances values from a search it should save a copy of them before running subsearches against the same set of points. If this is not done then the distance values from the first search may be overwritten by those of the subsequent searches.
Geo::Index uses the spherical haversine formula to compute distances. While quite fast, its accuracy over long distances is poor, with a worst case error of about 0.1% (22 km). Since the module already has provision for changing the backends used for the distance methods, adding a new backend to, for example, compute accurate distances on e.g. a WGS-84 spheroid would be simple and straight-forward.
In places the code can be repetitious or awkward in style. This was done because, especially in the inner loops, speed has been favoured over clarity.
Please submit any bugs or feature requests either to bug-geo-index at rt.cpan.org, through CPAN's web interface, or through Github. In any case I will receive notification when you do and you will be automatically notified of progress on your submission as it takes place. Any other comments can be sent to akh@cpan.org.
bug-geo-index at rt.cpan.org
akh@cpan.org
0.0.8 (2019-04-10) - Corrected internal links, CPAN release
Corrected POD links to use spaces instead of underscores
Uploaded to CPAN
0.0.7 (2019-04-08) - Methods can now be called using snake case
Added method aliases so that CamelCase methods can also be called using snake_case.
0.0.6 (2019-04-05) - Bug fixes, additional tests
GetIntLatLon(...): Fixed off-by-one error at the north pole
GetIntLatLon(...)
This affected Index(...) and Search(...).
GetIntLat(...): Fixed off-by-one error at the north pole
GetIntLat(...)
More thorough tests
Improved documentation
0.0.5 (2019-04-04) - Added methods, enhancements
PointCount( ): New method
PointCount( )
AllPoints( ): New method
AllPoints( )
Sweep(...): Added extra_keys option
extra_keys
Vacuum(...): Added extra_keys option
0.0.4 (2019-04-03) - Switched from Inline::C to XS for low-level C functions, minor restructuring
Low-level C code is now in Index.xs.
Index.xs
All references to Inline::C have been removed.
Deprecated DeletePointIndex(...) and replaced it with Unindex(...)
0.0.3 (2019-04-01) - Added Vacuum(...), Sweep(...), and tests plus bug fixes and minor enhancements
Sweep(...): New method
Vacuum(...): New method
Added tests
SearchByBounds(...): Bug fixes
new(...): Added quiet option
quiet
0.0.2 (2019-03-31) - Bug fixes and minor enhancements
Index(...): Fixed bug for points added near (but not at) the north pole
GetConfiguration( ): Added supported_key_types, supported_code_types, and tile_meters values>
0.0.1 (2018-03-30) - Initial release
Alex Kent Hajnal akh@cpan.org https://alephnull.net/software
Geo::Index
Copyright 2019 Alexander Hajnal, All rights reserved
This module is free software; you can redistribute it and/or modify it under the same terms as Perl itself. See perlartistic.
To install Geo::Index, copy and paste the appropriate command in to your terminal.
cpanm
cpanm Geo::Index
CPAN shell
perl -MCPAN -e shell install Geo::Index
For more information on module installation, please visit the detailed CPAN module installation guide.