Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Remove lonlat_2d and cartesian_2d types #662

Merged
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
15 changes: 7 additions & 8 deletions cpp/benchmarks/pairwise_linestring_distance.cu
Original file line number Diff line number Diff line change
Expand Up @@ -64,24 +64,23 @@ using namespace cuspatial;
*
*/
template <typename T>
std::tuple<rmm::device_vector<cartesian_2d<T>>, rmm::device_vector<int32_t>> generate_linestring(
int32_t num_strings, int32_t num_segments_per_string, T segment_length, cartesian_2d<T> init_xy)
std::tuple<rmm::device_vector<vec_2d<T>>, rmm::device_vector<int32_t>> generate_linestring(
int32_t num_strings, int32_t num_segments_per_string, T segment_length, vec_2d<T> init_xy)
{
int32_t num_points = num_strings * (num_segments_per_string + 1);

auto offset_iter = detail::make_counting_transform_iterator(
0, [num_segments_per_string](auto i) { return i * num_segments_per_string; });
auto rads_iter = detail::make_counting_transform_iterator(0, [](auto i) {
return cartesian_2d<T>{cos(static_cast<T>(i)), sin(static_cast<T>(i))};
return vec_2d<T>{cos(static_cast<T>(i)), sin(static_cast<T>(i))};
});

std::vector<int32_t> offsets(offset_iter, offset_iter + num_strings);
std::vector<cartesian_2d<T>> rads(rads_iter, rads_iter + num_points);
std::vector<cartesian_2d<T>> points(num_points);
std::vector<vec_2d<T>> rads(rads_iter, rads_iter + num_points);
std::vector<vec_2d<T>> points(num_points);

auto random_walk_func = [segment_length](cartesian_2d<T> const& prev,
cartesian_2d<T> const& rad) {
return cartesian_2d<T>{prev.x + segment_length * rad.x, prev.y + segment_length * rad.y};
auto random_walk_func = [segment_length](vec_2d<T> const& prev, vec_2d<T> const& rad) {
return prev + segment_length * rad;
};

thrust::exclusive_scan(
Expand Down
34 changes: 16 additions & 18 deletions cpp/benchmarks/point_in_polygon.cu
Original file line number Diff line number Diff line change
Expand Up @@ -39,11 +39,11 @@ auto constexpr num_rings_per_polygon = 1; // only 1 ring for now
* @brief Generate a random point within a window of [minXY, maxXY]
*/
template <typename T>
cartesian_2d<T> random_point(cartesian_2d<T> minXY, cartesian_2d<T> maxXY)
vec_2d<T> random_point(vec_2d<T> minXY, vec_2d<T> maxXY)
{
auto x = minXY.x + (maxXY.x - minXY.x) * rand() / static_cast<T>(RAND_MAX);
auto y = minXY.y + (maxXY.y - minXY.y) * rand() / static_cast<T>(RAND_MAX);
return cartesian_2d<T>{x, y};
return vec_2d<T>{x, y};
}

/**
Expand All @@ -63,14 +63,12 @@ cartesian_2d<T> random_point(cartesian_2d<T> minXY, cartesian_2d<T> maxXY)
*
*/
template <typename T>
std::tuple<rmm::device_vector<int32_t>,
rmm::device_vector<int32_t>,
rmm::device_vector<cartesian_2d<T>>>
generate_polygon(int32_t num_sides, T radius, cartesian_2d<T> minXY, cartesian_2d<T> maxXY)
std::tuple<rmm::device_vector<int32_t>, rmm::device_vector<int32_t>, rmm::device_vector<vec_2d<T>>>
generate_polygon(int32_t num_sides, T radius, vec_2d<T> minXY, vec_2d<T> maxXY)
{
std::vector<int32_t> polygon_offsets(num_polygons);
std::vector<int32_t> ring_offsets(num_polygons * num_rings_per_polygon);
std::vector<cartesian_2d<T>> polygon_points(31 * (num_sides + 1));
std::vector<vec_2d<T>> polygon_points(31 * (num_sides + 1));

std::iota(polygon_offsets.begin(), polygon_offsets.end(), 0);
std::iota(ring_offsets.begin(), ring_offsets.end(), 0);
Expand All @@ -84,11 +82,11 @@ generate_polygon(int32_t num_sides, T radius, cartesian_2d<T> minXY, cartesian_2
auto begin = i * num_sides + polygon_points.begin();
auto center = random_point(minXY, maxXY);
std::transform(it, it + num_sides + 1, begin, [num_sides, radius, center](int32_t j) {
return cartesian_2d<T>{
static_cast<T>(radius * std::cos(2 * PI * (j % num_sides) / static_cast<T>(num_sides)) +
center.x),
static_cast<T>(radius * std::sin(2 * PI * (j % num_sides) / static_cast<T>(num_sides)) +
center.y)};
return center +
radius *
vec_2d<T>{
static_cast<T>(std::cos(2 * PI * (j % num_sides) / static_cast<T>(num_sides))),
static_cast<T>(std::sin(2 * PI * (j % num_sides) / static_cast<T>(num_sides)))};
});
}

Expand All @@ -102,11 +100,11 @@ generate_polygon(int32_t num_sides, T radius, cartesian_2d<T> minXY, cartesian_2
* @tparam T The floating point type for the coordinates
*/
template <typename T>
rmm::device_vector<cartesian_2d<T>> generate_points(int32_t num_test_points,
cartesian_2d<T> minXY,
cartesian_2d<T> maxXY)
rmm::device_vector<vec_2d<T>> generate_points(int32_t num_test_points,
vec_2d<T> minXY,
vec_2d<T> maxXY)
{
std::vector<cartesian_2d<T>> points(num_test_points);
std::vector<vec_2d<T>> points(num_test_points);
std::generate(
points.begin(), points.end(), [minXY, maxXY]() { return random_point(minXY, maxXY); });
// Implicitly convert to device_vector
Expand All @@ -120,8 +118,8 @@ void point_in_polygon_benchmark(nvbench::state& state, nvbench::type_list<T>)
cuspatial::rmm_pool_raii rmm_pool;

std::srand(0); // For reproducibility
auto const minXY = cartesian_2d<T>{-radius * 2, -radius * 2};
auto const maxXY = cartesian_2d<T>{radius * 2, radius * 2};
auto const minXY = vec_2d<T>{-radius * 2, -radius * 2};
auto const maxXY = vec_2d<T>{radius * 2, radius * 2};

auto const num_test_points{state.get_int64("NumTestPoints")},
num_sides_per_ring{state.get_int64("NumSidesPerRing")};
Expand Down
56 changes: 27 additions & 29 deletions cpp/doc/libcuspatial_refactoring_guide.md
Original file line number Diff line number Diff line change
Expand Up @@ -48,8 +48,8 @@ There are a few key points to notice.
1. The API is very similar to STL algorithms such as `std::transform`.
2. All array inputs and outputs are iterator type templates.
3. Longitude/Latitude data is passed as array of structures, using the `cuspatial::lonlat_2d`
type (include/cuspatial/types.hpp). This is enforced using a `static_assert` in the function
3. Longitude/Latitude data is passed as array of structures, using the `cuspatial::vec_2d`
type (include/cuspatial/vec_2d.hpp). This is enforced using a `static_assert` in the function
body (discussed later).
4. The `Location` type is a template that is by default equal to the `value_type` of the input
iterators.
Expand Down Expand Up @@ -96,8 +96,7 @@ Following is the (Doxygen) documentation for the above `cuspatial::haversine_dis
* [LegacyRandomAccessIterator][LinkLRAI] and be device-accessible.
* @tparam OutputIt Output iterator. Must meet the requirements of
* [LegacyRandomAccessIterator][LinkLRAI] and be device-accessible.
* @tparam Location The `value_type` of `LonLatItA` and `LonLatItB`. Must be
* `cuspatial::lonlat_2d<T>`.
* @tparam Location The `value_type` of `LonLatItA` and `LonLatItB`. Must be `cuspatial::vec_2d<T>`.
* @tparam T The underlying coordinate type. Must be a floating-point type.
*
* @pre `a_lonlat_first` may equal `distance_first`, but the range `[a_lonlat_first, a_lonlat_last)`
Expand All @@ -107,7 +106,7 @@ Following is the (Doxygen) documentation for the above `cuspatial::haversine_dis
* shall not overlap the range `[distance_first, distance_first + (b_lonlat_last - b_lonlat_last))
* otherwise.
* @pre All iterators must have the same `Location` type, with the same underlying floating-point
* coordinate type (e.g. `cuspatial::lonlat_2d<float>`).
* coordinate type (e.g. `cuspatial::vec_2d<float>`).
*
* @return Output iterator to the element past the last distance computed.
*
Expand All @@ -132,13 +131,16 @@ This is the existing API, unchanged by refactoring. Here is the existing
`cuspatial::haversine_distance`:

```C++
std::unique_ptr<cudf::column> haversine_distance(
cudf::column_view const& a_lon,
cudf::column_view const& a_lat,
cudf::column_view const& b_lon,
cudf::column_view const& b_lat,
double const radius = EARTH_RADIUS_KM,
rmm::mr::device_memory_resource* mr = rmm::mr::get_current_device_resource());
template <class LonLatItA,
class LonLatItB,
class OutputIt,
class T = typename detail::iterator_vec_base_type<LonLatItA>>
OutputIt haversine_distance(LonLatItA a_lonlat_first,
LonLatItA a_lonlat_last,
LonLatItB b_lonlat_first,
OutputIt distance_first,
T const radius = EARTH_RADIUS_KM,
rmm::cuda_stream_view stream = rmm::cuda_stream_default);
```
key points:
Expand Down Expand Up @@ -186,27 +188,23 @@ a simple matter of a few static asserts and dynamic expectation checks, followed
`thrust::transform` with a custom transform functor.

```C++
template <class LonLatItA,
class LonLatItB,
class OutputIt,
class Location = typename std::iterator_traits<LonLatItA>::value_type,
class T = typename Location::value_type>
template <class LonLatItA, class LonLatItB, class OutputIt, class T>
OutputIt haversine_distance(LonLatItA a_lonlat_first,
LonLatItA a_lonlat_last,
LonLatItB b_lonlat_first,
OutputIt distance_first,
T const radius,
rmm::cuda_stream_view stream)
{
using LocationB = typename std::iterator_traits<LonLatItB>::value_type;
static_assert(std::conjunction_v<std::is_same<lonlat_2d<T>, Location>,
std::is_same<lonlat_2d<T>, LocationB>>,
"Inputs must be cuspatial::lonlat_2d");
static_assert(detail::is_same<vec_2d<T>,
detail::iterator_value_type<LonLatItA>,
detail::iterator_value_type<LonLatItB>>(),
"Inputs must be cuspatial::vec_2d");
static_assert(
std::conjunction_v<std::is_floating_point<T>,
std::is_floating_point<typename LocationB::value_type>,
std::is_floating_point<typename std::iterator_traits<OutputIt>::value_type>>,
"Haversine distance supports only floating-point coordinates.");
detail::is_same_floating_point<T,
typename detail::iterator_vec_base_type<LonLatItA>,
typename detail::iterator_value_type<OutputIt>>(),
"All iterator types and radius must have the same floating-point coordinate value type.");

CUSPATIAL_EXPECTS(radius > 0, "radius must be positive.");

Expand All @@ -221,15 +219,15 @@ OutputIt haversine_distance(LonLatItA a_lonlat_first,
Note that we `static_assert` that the types of the iterator inputs match documented expectations.
We also do a runtime check that the radius is positive. Finally we just call `thrust::transform`,
passing it an instance of `haversine_distance_functor`, which is a function of two `lonlat_2d`
passing it an instance of `haversine_distance_functor`, which is a function of two `vec_2d<T>`
inputs that implements the Haversine distance formula.
### libcudf-based API Implementation
The substance of the refactoring is making the libcudf-based API a wrapper around the header-only
API. This mostly involves replacing business logic implementation in the type-dispatched functor
with a call to the header-only API. We also need to convert disjoint latitude and longitude inputs
into `lonlat_2d<T>` structs. This is easily done using the `cuspatial::make_lonlat_iterator` utility
into `vec_2d<T>` structs. This is easily done using the `cuspatial::make_vec_2d_iterator` utility
provided in `type_utils.hpp`.
So, to refactor the libcudf-based API, we remove the following code.
Expand Down Expand Up @@ -259,8 +257,8 @@ thrust::transform(rmm::exec_policy(stream),
And replace it with the following code.

```C++
auto lonlat_a = cuspatial::make_lonlat_iterator(a_lon.begin<T>(), a_lat.begin<T>());
auto lonlat_b = cuspatial::make_lonlat_iterator(b_lon.begin<T>(), b_lat.begin<T>());
auto lonlat_a = cuspatial::make_vec_2d_iterator(a_lon.begin<T>(), a_lat.begin<T>());
auto lonlat_b = cuspatial::make_vec_2d_iterator(b_lon.begin<T>(), b_lat.begin<T>());

cuspatial::haversine_distance(lonlat_a,
lonlat_a + a_lon.size(),
Expand Down
1 change: 1 addition & 0 deletions cpp/include/cuspatial/detail/utility/traits.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -16,6 +16,7 @@

#pragma once

#include <iterator>
#include <type_traits>

namespace cuspatial {
Expand Down
8 changes: 4 additions & 4 deletions cpp/include/cuspatial/experimental/coordinate_transform.cuh
Original file line number Diff line number Diff line change
Expand Up @@ -33,9 +33,9 @@ namespace cuspatial {
* @param[out] xy_first: beginning of range of output x/y coordinates.
* @param[in] stream: The CUDA stream on which to perform computations and allocate memory.
*
* All input iterators must have a `value_type` of `cuspatial::lonlat_2d<T>` (Lat/Lon coordinates),
* and the output iterator must be able to accept for storage values of type
* `cuspatial::cartesian_2d<T>` (Cartesian coordinates).
* All input iterators must have a `value_type` of `cuspatial::vec_2d<T>` (Lat/Lon coordinates),
* and the output iterator must be able to accept for storage values of type `cuspatial::vec_2d<T>`
* (Cartesian coordinates).
*
* @tparam InputIt Iterator over longitude/latitude locations. Must meet the requirements of
* [LegacyRandomAccessIterator][LinkLRAI] and be device-accessible.
Expand All @@ -56,7 +56,7 @@ template <class InputIt, class OutputIt, class T>
OutputIt lonlat_to_cartesian(InputIt lon_lat_first,
InputIt lon_lat_last,
OutputIt xy_first,
lonlat_2d<T> origin,
vec_2d<T> origin,
rmm::cuda_stream_view stream = rmm::cuda_stream_default);

} // namespace cuspatial
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -53,16 +53,16 @@ __device__ inline T lat_to_y(T lat)

template <typename T>
struct to_cartesian_functor {
to_cartesian_functor(lonlat_2d<T> origin) : _origin(origin) {}
to_cartesian_functor(vec_2d<T> origin) : _origin(origin) {}

cartesian_2d<T> __device__ operator()(lonlat_2d<T> loc)
vec_2d<T> __device__ operator()(vec_2d<T> loc)
{
return cartesian_2d<T>{lon_to_x(_origin.x - loc.x, midpoint(loc.y, _origin.y)),
lat_to_y(_origin.y - loc.y)};
return vec_2d<T>{lon_to_x(_origin.x - loc.x, midpoint(loc.y, _origin.y)),
lat_to_y(_origin.y - loc.y)};
}

private:
lonlat_2d<T> _origin{};
vec_2d<T> _origin{};
};

} // namespace detail
Expand All @@ -71,7 +71,7 @@ template <class InputIt, class OutputIt, class T>
OutputIt lonlat_to_cartesian(InputIt lon_lat_first,
InputIt lon_lat_last,
OutputIt xy_first,
lonlat_2d<T> origin,
vec_2d<T> origin,
rmm::cuda_stream_view stream)
{
static_assert(std::is_floating_point_v<T>,
Expand Down
30 changes: 15 additions & 15 deletions cpp/include/cuspatial/experimental/detail/haversine.cuh
Original file line number Diff line number Diff line change
Expand Up @@ -17,14 +17,14 @@
#pragma once

#include <cuspatial/constants.hpp>
#include <cuspatial/detail/utility/traits.hpp>
#include <cuspatial/error.hpp>
#include <cuspatial/vec_2d.hpp>

#include <rmm/cuda_stream_view.hpp>
#include <rmm/exec_policy.hpp>

#include <thrust/transform.h>
#include <thrust/tuple.h>

#include <type_traits>

Expand All @@ -36,12 +36,12 @@ template <typename T>
struct haversine_distance_functor {
haversine_distance_functor(T radius) : radius_(radius) {}

__device__ T operator()(lonlat_2d<T> point_a, lonlat_2d<T> point_b)
__device__ T operator()(vec_2d<T> lonlat_a, vec_2d<T> lonlat_b)
{
auto ax = point_a.x * DEGREE_TO_RADIAN;
auto ay = point_a.y * DEGREE_TO_RADIAN;
auto bx = point_b.x * DEGREE_TO_RADIAN;
auto by = point_b.y * DEGREE_TO_RADIAN;
auto ax = lonlat_a.x * DEGREE_TO_RADIAN;
auto ay = lonlat_a.y * DEGREE_TO_RADIAN;
auto bx = lonlat_b.x * DEGREE_TO_RADIAN;
auto by = lonlat_b.y * DEGREE_TO_RADIAN;

// haversine formula
auto x = (bx - ax) / 2;
Expand All @@ -58,23 +58,23 @@ struct haversine_distance_functor {

} // namespace detail

template <class LonLatItA, class LonLatItB, class OutputIt, class Location, class T>
template <class LonLatItA, class LonLatItB, class OutputIt, class T>
OutputIt haversine_distance(LonLatItA a_lonlat_first,
LonLatItA a_lonlat_last,
LonLatItB b_lonlat_first,
OutputIt distance_first,
T const radius,
rmm::cuda_stream_view stream)
{
using LocationB = typename std::iterator_traits<LonLatItB>::value_type;
static_assert(detail::is_same<vec_2d<T>,
detail::iterator_value_type<LonLatItA>,
detail::iterator_value_type<LonLatItB>>(),
"Inputs must be cuspatial::vec_2d");
static_assert(
std::conjunction_v<std::is_same<lonlat_2d<T>, Location>, std::is_same<lonlat_2d<T>, LocationB>>,
"Inputs must be cuspatial::lonlat_2d");
static_assert(
std::conjunction_v<std::is_floating_point<T>,
std::is_floating_point<typename LocationB::value_type>,
std::is_floating_point<typename std::iterator_traits<OutputIt>::value_type>>,
"Haversine distance supports only floating-point coordinates.");
detail::is_same_floating_point<T,
typename detail::iterator_vec_base_type<LonLatItA>,
typename detail::iterator_value_type<OutputIt>>(),
"All iterator types and radius must have the same floating-point coordinate value type.");

CUSPATIAL_EXPECTS(radius > 0, "radius must be positive.");

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -148,21 +148,17 @@ void pairwise_linestring_distance(OffsetIterator linestring1_offsets_first,
OutputIt distances_first,
rmm::cuda_stream_view stream)
{
using T = typename std::iterator_traits<Cart2dItA>::value_type::value_type;

static_assert(
detail::is_floating_point<T,
typename std::iterator_traits<Cart2dItB>::value_type::value_type,
typename std::iterator_traits<OutputIt>::value_type>(),
"Inputs and output must have floating point value type.");
using T = typename detail::iterator_vec_base_type<Cart2dItA>;

static_assert(detail::is_same<T, typename std::iterator_traits<OutputIt>::value_type>(),
"Inputs and output must have the same value type.");
static_assert(detail::is_same_floating_point<T,
typename detail::iterator_vec_base_type<Cart2dItB>,
typename detail::iterator_value_type<OutputIt>>(),
"Inputs and output must have the same floating point value type.");

static_assert(detail::is_same<cartesian_2d<T>,
typename std::iterator_traits<Cart2dItA>::value_type,
typename std::iterator_traits<Cart2dItB>::value_type>(),
"All input types must be cuspatial::cartesian_2d with the same value type");
static_assert(detail::is_same<vec_2d<T>,
typename detail::iterator_value_type<Cart2dItA>,
typename detail::iterator_value_type<Cart2dItB>>(),
"All input types must be cuspatial::vec_2d with the same value type");

auto const num_linestring_pairs =
thrust::distance(linestring1_offsets_first, linestring1_offsets_last);
Expand Down
Loading