Trajectory#

Functions for identifying and grouping trajectories from point data.

cuspatial.derive_trajectories(object_ids, xs, ys, timestamps)#

Derive trajectories from object ids, points, and timestamps.

Parameters
object_ids

column of object (e.g., vehicle) ids

xs

column of x-coordinates (in kilometers)

ys

column of y-coordinates (in kilometers)

timestamps

column of timestamps in any resolution

Returns
resulttuple (objects, traj_offsets)
objectscudf.DataFrame

object_ids, xs, ys, and timestamps sorted by (object_id, timestamp), used by trajectory_bounding_boxes and trajectory_distances_and_speeds

traj_offsetscudf.Series

offsets of discovered trajectories

Examples

Compute sorted objects and discovered trajectories

>>> objects, traj_offsets = cuspatial.derive_trajectories(
        [0, 1, 2, 3],  # object_id
        [0, 0, 1, 1],  # x
        [0, 0, 1, 1],  # y
        [0, 10, 0, 10] # timestamp
    )
>>> print(traj_offsets)
    0  0
    1  2
>>> print(objects)
       object_id       x       y  timestamp
    0          0       1       0          0
    1          0       0       0         10
    2          1       3       1          0
    3          1       2       1         10
cuspatial.trajectory_distances_and_speeds(num_trajectories, object_ids, xs, ys, timestamps)#

Compute the distance traveled and speed of sets of trajectories

Parameters
num_trajectories

number of trajectories (unique object ids)

object_ids

column of object (e.g., vehicle) ids

xs

column of x-coordinates (in kilometers)

ys

column of y-coordinates (in kilometers)

timestamps

column of timestamps in any resolution

Returns
resultcudf.DataFrame
meterscudf.Series

trajectory distance (in kilometers)

speedcudf.Series

trajectory speed (in meters/second)

Examples

Compute the distances and speeds of derived trajectories

>>> objects, traj_offsets = cuspatial.derive_trajectories(...)
>>> dists_and_speeds = cuspatial.trajectory_distances_and_speeds(
        len(traj_offsets)
        objects['object_id'],
        objects['x'],
        objects['y'],
        objects['timestamp']
    )
>>> print(dists_and_speeds)
                   distance          speed
    trajectory_id
    0                1000.0  100000.000000
    1                1000.0  111111.109375
cuspatial.directed_hausdorff_distance(xs, ys, space_offsets)#

Compute the directed Hausdorff distances between all pairs of spaces.

Parameters
xs

column of x-coordinates

ys

column of y-coordinates

space_offsets

beginning index of each space, plus the last space’s end offset.

Returns
resultcudf.DataFrame

The pairwise directed distance matrix with one row and one column per input space; the value at row i, column j represents the hausdorff distance from space i to space j.

Examples

The directed Hausdorff distance from one space to another is the greatest of all the distances between any point in the first space to the closest point in the second.

Wikipedia

Consider a pair of lines on a grid:

     :
     x
-----xyy---
     :
     :

x0 = (0, 0), x1 = (0, 1)

y0 = (1, 0), y1 = (2, 0)

x0 is the closest point in x to y. The distance from x0 to the farthest point in y is 2.

y0 is the closest point in y to x. The distance from y0 to the farthest point in x is 1.414.

Compute the directed hausdorff distances between a set of spaces

>>> result = cuspatial.directed_hausdorff_distance(
        [0, 1, 0, 0], # xs
        [0, 0, 1, 2], # ys
        [0, 2, 4],    # space_offsets
    )
>>> print(result)
         0         1
    0  0.0  1.414214
    1  2.0  0.000000
cuspatial.trajectory_bounding_boxes(num_trajectories, object_ids, xs, ys)#

Compute the bounding boxes of sets of trajectories.

Parameters
num_trajectories

number of trajectories (unique object ids)

object_ids

column of object (e.g., vehicle) ids

xs

column of x-coordinates (in kilometers)

ys

column of y-coordinates (in kilometers)

Returns
resultcudf.DataFrame

minimum bounding boxes (in kilometers) for each trajectory

x_mincudf.Series

the minimum x-coordinate of each bounding box

y_mincudf.Series

the minimum y-coordinate of each bounding box

x_maxcudf.Series

the maximum x-coordinate of each bounding box

y_maxcudf.Series

the maximum y-coordinate of each bounding box

Examples

Compute the minimum bounding boxes of derived trajectories

>>> objects, traj_offsets = trajectory.derive_trajectories(
        [0, 0, 1, 1],  # object_id
        [0, 1, 2, 3],  # x
        [0, 0, 1, 1],  # y
        [0, 10, 0, 10] # timestamp
    )
>>> traj_bounding_boxes = cuspatial.trajectory_bounding_boxes(
        len(traj_offsets),
        objects['object_id'],
        objects['x'],
        objects['y']
    )
>>> print(traj_bounding_boxes)
    x_min   y_min   x_max   y_max
0     0.0     0.0     2.0     2.0
1     1.0     1.0     3.0     3.0
class cuspatial.CubicSpline(t, y, ids=None, size=None, prefixes=None)#

Fits each column of the input Series y to a hermetic cubic spline.

cuspatial.CubicSpline supports two usage patterns: The first is identical to scipy.interpolate.CubicSpline:

curve = cuspatial.CubicSpline(t, y)
new_points = curve(np.linspace(t.min, t.max, 50))

This allows API parity with scipy. This isn’t recommended, as scipy host based interpolation performance is likely to exceed GPU performance for a single curve.

However, cuSpatial massively outperforms scipy when many splines are fit simultaneously. Data must be arranged in a SoA format, and the exclusive prefix_sum of the separate curves must also be passed to the function.:

NUM_SPLINES = 100000
SPLINE_LENGTH = 101
t = cudf.Series(
    np.hstack((np.arange(SPLINE_LENGTH),) * NUM_SPLINES)
).astype('float32')
y = cudf.Series(
    np.random.random(SPLINE_LENGTH*NUM_SPLINES)
).astype('float32')
prefix_sum = cudf.Series(
    cp.arange(NUM_SPLINES + 1)*SPLINE_LENGTH
).astype('int32')
curve = cuspatial.CubicSpline(t, y, prefixes=prefix_sum)
new_samples = cudf.Series(
    np.hstack((np.linspace(
        0, (SPLINE_LENGTH - 1), (SPLINE_LENGTH - 1) * 2 + 1
    ),) * NUM_SPLINES)
).astype('float32')
curve_ids = cudf.Series(np.repeat(
    np.arange(0, NUM_SPLINES), SPLINE_LENGTH * 2 - 1
), dtype="int32")
new_points = curve(new_samples, curve_ids)

Methods

__call__(coordinates[, groups])

Interpolates new input values coordinates using the .c DataFrame or map of DataFrames.

CubicSpline.__init__(t, y, ids=None, size=None, prefixes=None)#

Computes various error preconditions on the input data, then uses CUDA to compute cubic splines for each set of input coordinates on the GPU in parallel.

Parameters
tcudf.Series

time sample values. Must be monotonically increasing.

ycudf.Series

columns to have curves fit to according to x

ids (Optional)cudf.Series

ids of each spline

size (Optional)cudf.Series

fixed size of each spline

prefixes (Optional)cudf.Series

alternative to size, allows splines of varying length. Not yet fully supported.

Returns
CubicSplinecallable o

o.c contains the coefficients that can be used to compute new points along the spline fitting the original t data. o(n) interpolates the spline coordinates along new input values n.

CubicSpline.__call__(coordinates, groups=None)#

Interpolates new input values coordinates using the .c DataFrame or map of DataFrames.