Remove lonlat_2d and cartesian_2d types (#662)

Closes #653.  Removes `lonlat_2d<T>` and `cartesian_2d` structs and related iterator factories and replaces all usage of these structs with their base class `vec_2d<T>. See #653 for discussion / motivation.

Authors:
  - Mark Harris (https://github.com/harrism)

Approvers:
  - H. Thomson Comer (https://github.com/thomcom)
  - Michael Wang (https://github.com/isVoid)
  - Paul Taylor (https://github.com/trxcllnt)

URL: #662

cpp/benchmarks/pairwise_linestring_distance.cu

@@ -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(

cpp/benchmarks/point_in_polygon.cu

@@ -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};
 }
 
 /**
@@ -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);
@@ -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)))};
     });
   }
 
@@ -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
@@ -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")};

cpp/doc/libcuspatial_refactoring_guide.md

@@ -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.
@@ -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)`
@@ -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.
  *
@@ -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:
@@ -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.");
 
@@ -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.
@@ -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(),

cpp/include/cuspatial/detail/utility/traits.hpp

@@ -16,6 +16,7 @@
 
 #pragma once
 
+#include <iterator>
 #include <type_traits>
 
 namespace cuspatial {

cpp/include/cuspatial/experimental/coordinate_transform.cuh

@@ -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.
@@ -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

cpp/include/cuspatial/experimental/detail/coordinate_transform.cuh

@@ -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
@@ -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>,

cpp/include/cuspatial/experimental/detail/haversine.cuh

@@ -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>
 
@@ -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;
@@ -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.");
 

cpp/include/cuspatial/experimental/detail/linestring_distance.cuh

@@ -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);