Vectorized Load, refactors type_utils.hpp into iterator_factory.cuh

cpp/benchmarks/pairwise_linestring_distance.cu

@@ -18,8 +18,8 @@
 #include <nvbench/nvbench.cuh>
 
 #include <cuspatial/detail/iterator.hpp>
+#include <cuspatial/experimental/iterator_factory.cuh>
 #include <cuspatial/experimental/linestring_distance.cuh>
-#include <cuspatial/experimental/type_utils.hpp>
 #include <cuspatial/vec_2d.hpp>
 
 #include <rmm/device_vector.hpp>

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

@@ -88,17 +88,18 @@ template <class Cart2dItA,
           class OffsetIteratorB,
           class OffsetIteratorC,
           class OutputIterator>
-void __global__ pairwise_point_linestring_distance(OffsetIteratorA point_geometry_offset_first,
-                                                   OffsetIteratorA point_geometry_offset_last,
-                                                   Cart2dItA points_first,
-                                                   Cart2dItA points_last,
-                                                   OffsetIteratorB linestring_geometry_offset_first,
-                                                   OffsetIteratorB linestring_geometry_offset_last,
-                                                   OffsetIteratorC linestring_part_offsets_first,
-                                                   OffsetIteratorC linestring_part_offsets_last,
-                                                   Cart2dItB linestring_points_first,
-                                                   Cart2dItB linestring_points_last,
-                                                   OutputIterator distances)
+void __global__
+pairwise_point_linestring_distance_kernel(OffsetIteratorA point_geometry_offset_first,
+                                          OffsetIteratorA point_geometry_offset_last,
+                                          Cart2dItA points_first,
+                                          Cart2dItA points_last,
+                                          OffsetIteratorB linestring_geometry_offset_first,
+                                          OffsetIteratorB linestring_geometry_offset_last,
+                                          OffsetIteratorC linestring_part_offsets_first,
+                                          OffsetIteratorC linestring_part_offsets_last,
+                                          Cart2dItB linestring_points_first,
+                                          Cart2dItB linestring_points_last,
+                                          OutputIterator distances)
 {
   using T = iterator_vec_base_type<Cart2dItA>;
 
@@ -193,18 +194,19 @@ OutputIt pairwise_point_linestring_distance(OffsetIteratorA point_geometry_offse
   std::size_t const num_blocks =
     (num_linestring_points + threads_per_block - 1) / threads_per_block;
 
-  detail::pairwise_point_linestring_distance<<<num_blocks, threads_per_block, 0, stream.value()>>>(
-    point_geometry_offset_first,
-    point_geometry_offset_last,
-    points_first,
-    points_last,
-    linestring_geometry_offset_first,
-    linestring_geometry_offset_first + num_pairs + 1,
-    linestring_part_offsets_first,
-    linestring_part_offsets_last,
-    linestring_points_first,
-    linestring_points_last,
-    distances_first);
+  detail::
+    pairwise_point_linestring_distance_kernel<<<num_blocks, threads_per_block, 0, stream.value()>>>(
+      point_geometry_offset_first,
+      point_geometry_offset_last,
+      points_first,
+      points_last,
+      linestring_geometry_offset_first,
+      linestring_geometry_offset_first + num_pairs + 1,
+      linestring_part_offsets_first,
+      linestring_part_offsets_last,
+      linestring_points_first,
+      linestring_points_last,
+      distances_first);
 
   CUSPATIAL_CUDA_TRY(cudaGetLastError());
 

cpp/include/cuspatial/experimental/type_utils.hpp → cpp/include/cuspatial/experimental/iterator_factory.cuh

@@ -14,22 +14,43 @@
  * limitations under the License.
  */
 
+#pragma once
+
 #include <cuspatial/detail/iterator.hpp>
 #include <cuspatial/traits.hpp>
 #include <cuspatial/vec_2d.hpp>
 
+#include <thrust/detail/raw_reference_cast.h>
 #include <thrust/iterator/permutation_iterator.h>
 #include <thrust/iterator/transform_iterator.h>
 #include <thrust/iterator/transform_output_iterator.h>
 #include <thrust/iterator/zip_iterator.h>
 #include <thrust/tuple.h>
+#include <thrust/type_traits/is_contiguous_iterator.h>
 
 #include <type_traits>
 
 namespace cuspatial {
-
 namespace detail {
 
+/**
+ * @brief Get corresponding 2-element vector type from `Element` type.
+ * Only for floating point types.
+ */
+template <typename Element>
+struct element_to_element2 {
+};
+
+template <>
+struct element_to_element2<float> {
+  using type = float2;
+};
+
+template <>
+struct element_to_element2<double> {
+  using type = double2;
+};
+
 /**
  * @internal
  * @brief Helper to convert a tuple of elements into a `vec_2d`
@@ -54,6 +75,74 @@ struct vec_2d_to_tuple {
   }
 };
 
+/**
+ * @internal
+ * @brief Generic to convert any iterator pointing to interleaved xy range into
+ * iterator of vec_2d.
+ *
+ * This generic version does not make use of vectorized load.
+ *
+ * @pre `Iter` has operator[] defined.
+ * @pre std::iterator_traits<Iter>::value_type is convertible to `T`.
+ */
+template <typename Iter, typename Enable = void>
+struct interleaved_to_vec_2d {
+  using T = typename std::iterator_traits<Iter>::value_type;
+  Iter it;
+  __device__ vec_2d<T> operator()(std::size_t i) { return vec_2d<T>{it[2 * i], it[2 * i + 1]}; }
+};
+
+/**
+ * @brief Specialization for raw pointers to interleaved xy range.
+ *
+ * @pre `Iter` is raw pointer
+ */
+template <typename Iter>
+struct interleaved_to_vec_2d<Iter, typename std::enable_if_t<std::is_pointer_v<Iter>>> {
+  using pointer_t = typename std::pointer_traits<Iter>::pointer;
+  using T         = typename std::remove_const_t<typename std::pointer_traits<Iter>::element_type>;
+  using T2        = typename element_to_element2<T>::type;
+  pointer_t ptr;
+
+  __device__ vec_2d<T> operator()(std::size_t i)
+  {
+    T2 const* f2it = reinterpret_cast<T2 const*>(ptr);
+    return vec_2d<T>{f2it[i].x, f2it[i].y};
+  }
+};
+
+/**
+ * @brief Specialization for thrust iterators conforming to `contiguous_iterator`.
+ *
+ * @pre `Iter` is a thrust iterator and conforms to `contiguous_iterator`.
+ */
+template <typename Iter>
+struct interleaved_to_vec_2d<
+  Iter,
+  typename std::enable_if_t<!std::is_pointer_v<Iter> && thrust::is_contiguous_iterator_v<Iter>>> {
+  using T  = typename std::iterator_traits<Iter>::value_type;
+  using T2 = typename element_to_element2<T>::type;
+  T* ptr;
+
+  interleaved_to_vec_2d(Iter it) { ptr = &thrust::raw_reference_cast(*it); }
+
+  __device__ vec_2d<T> operator()(std::size_t i)
+  {
+    T2 const* f2it = reinterpret_cast<T2 const*>(ptr);
+    return vec_2d<T>{f2it[i].x, f2it[i].y};
+  }
+};
+
+/**
+ * @internal
+ * @brief Functor to transform an index to strided index.
+ */
+struct strided_functor {
+  std::size_t _stride;
+  strided_functor(std::size_t stride) : _stride(stride) {}
+  auto __device__ operator()(std::size_t i) { return i * _stride; }
+};
+
 }  // namespace detail
 
 /**
@@ -64,14 +153,14 @@ struct vec_2d_to_tuple {
 /**
  * @brief Create an iterator to `vec_2d` data from two input iterators.
  *
- * Interleaves x and y coordinates from separate iterators into a single iterator to x-y
+ * Interleaves x and y coordinates from separate iterators into a single iterator to xy-
  * coordinates.
  *
  * @tparam VectorType cuSpatial vector type, must be `vec_2d`
- * @tparam FirstIter Iterator type to the first component of `vec_2d`. Must meet the requirements of
- * [LegacyRandomAccessIterator][LinkLRAI] and be device-accessible.
- * @tparam SecondIter Iterator type to the second component of `vec_2d`. Must meet the requirements
+ * @tparam FirstIter Iterator type to the first component of `vec_2d`. Must meet the requirements
  * of [LegacyRandomAccessIterator][LinkLRAI] and be device-accessible.
+ * @tparam SecondIter Iterator type to the second component of `vec_2d`. Must meet the
+ * requirements of [LegacyRandomAccessIterator][LinkLRAI] and be device-accessible.
  * @param first Iterator to beginning of `vec_2d::x`
  * @param second Iterator to beginning of `vec_2d::y`
  * @return Iterator to `vec_2d`
@@ -91,6 +180,21 @@ auto make_vec_2d_iterator(FirstIter first, SecondIter second)
   return thrust::make_transform_iterator(zipped, detail::tuple_to_vec_2d<T>());
 }
 
+/**
+ * @brief Create an iterator to `vec_2d` data from a single iterator.
+ *
+ * Creates a vec2d view from an interator to the starting range of interleaved x-y coordinates.
+ *
+ * @tparam
+ * @param d_points_begin
+ * @return
+ */
+template <typename Iter>
+auto make_vec_2d_iterator(Iter xy_begin)
+{
+  return detail::make_counting_transform_iterator(0, detail::interleaved_to_vec_2d<Iter>{xy_begin});
+}
+
 /**
  * @brief Create an output iterator to `vec_2d` data from two output iterators.
  *
@@ -120,25 +224,6 @@ auto make_zipped_vec_2d_output_iterator(FirstIter first, SecondIter second)
   return thrust::make_transform_output_iterator(zipped_out, detail::vec_2d_to_tuple<T>());
 }
 
-template <typename T, typename VectorType = vec_2d<T>>
-struct interleaved_to_vec_2d {
-  T const* it;
-  __device__ VectorType operator()(std::size_t i) { return VectorType{it[2 * i], it[2 * i + 1]}; }
-};
-
-template <typename Iter>
-auto interleaved_iterator_to_vec_2d_iterator(Iter d_points_begin)
-{
-  using T = typename std::iterator_traits<Iter>::value_type;
-  return detail::make_counting_transform_iterator(0, interleaved_to_vec_2d<T>{d_points_begin});
-}
-
-struct strided_functor {
-  std::size_t _stride;
-  strided_functor(std::size_t stride) : _stride(stride) {}
-  auto __device__ operator()(std::size_t i) { return i * _stride; }
-};
-
 /**
  * @brief Create an output iterator to `vec_2d` data from an iterator to an interleaved array.
  *
@@ -151,7 +236,7 @@ template <typename Iter>
 auto vec_2d_iterator_to_output_interleaved_iterator(Iter d_points_begin)
 {
   using T                     = typename std::iterator_traits<Iter>::value_type;
-  auto fixed_stride_2_functor = strided_functor(2);
+  auto fixed_stride_2_functor = detail::strided_functor(2);
   auto even_positions         = thrust::make_permutation_iterator(
     d_points_begin, detail::make_counting_transform_iterator(0, fixed_stride_2_functor));
   auto odd_positions = thrust::make_permutation_iterator(

cpp/include/doxygen_groups.h

@@ -108,7 +108,7 @@
  *          CuSpatial functions inside `experimental` folder are header-only and only accepts
  *          input/output iterators on coordinates. These factory functions are convenient ways
  *          to create iterators from data in various format.
- *          @file type_utils.hpp
+ *          @file iterator_factory.hpp
  *      @}
  * @}
  * @defgroup io I/O

cpp/src/spatial/hausdorff.cu

@@ -16,7 +16,7 @@
 
 #include <cuspatial/error.hpp>
 #include <cuspatial/experimental/hausdorff.cuh>
-#include <cuspatial/experimental/type_utils.hpp>
+#include <cuspatial/experimental/iterator_factory.cuh>
 
 #include <cudf/column/column.hpp>
 #include <cudf/column/column_device_view.cuh>

cpp/src/spatial/haversine.cu

@@ -17,7 +17,7 @@
 #include <cuspatial/constants.hpp>
 #include <cuspatial/error.hpp>
 #include <cuspatial/experimental/haversine.cuh>
-#include <cuspatial/experimental/type_utils.hpp>
+#include <cuspatial/experimental/iterator_factory.cuh>
 
 #include <cudf/column/column.hpp>
 #include <cudf/column/column_view.hpp>

cpp/src/spatial/linestring_distance.cu

@@ -15,8 +15,8 @@
  */
 
 #include <cuspatial/error.hpp>
+#include <cuspatial/experimental/iterator_factory.cuh>
 #include <cuspatial/experimental/linestring_distance.cuh>
-#include <cuspatial/experimental/type_utils.hpp>
 #include <cuspatial/vec_2d.hpp>
 
 #include <cudf/column/column_factories.hpp>

cpp/src/spatial/lonlat_to_cartesian.cu

@@ -16,7 +16,7 @@
 
 #include <cuspatial/error.hpp>
 #include <cuspatial/experimental/coordinate_transform.cuh>
-#include <cuspatial/experimental/type_utils.hpp>
+#include <cuspatial/experimental/iterator_factory.cuh>
 #include <cuspatial/vec_2d.hpp>
 
 #include <cudf/column/column.hpp>

cpp/src/spatial/point_distance.cu

@@ -15,8 +15,8 @@
  */
 
 #include <cuspatial/error.hpp>
+#include <cuspatial/experimental/iterator_factory.cuh>
 #include <cuspatial/experimental/point_distance.cuh>
-#include <cuspatial/experimental/type_utils.hpp>
 #include <cuspatial/vec_2d.hpp>
 
 #include <cudf/column/column_factories.hpp>
@@ -51,11 +51,11 @@ struct pairwise_point_distance_functor {
     auto points1_it = make_vec_2d_iterator(points1_x.begin<T>(), points1_y.begin<T>());
     auto points2_it = make_vec_2d_iterator(points2_x.begin<T>(), points2_y.begin<T>());
 
-    pairwise_point_distance(points1_it,
-                            points1_it + points1_x.size(),
-                            points2_it,
-                            distances->mutable_view().begin<T>(),
-                            stream);
+    cuspatial::pairwise_point_distance(points1_it,
+                                       points1_it + points1_x.size(),
+                                       points2_it,
+                                       distances->mutable_view().begin<T>(),
+                                       stream);
 
     return distances;
   }

cpp/src/spatial/point_in_polygon.cu

@@ -15,8 +15,8 @@
  */
 
 #include <cuspatial/error.hpp>
+#include <cuspatial/experimental/iterator_factory.cuh>
 #include <cuspatial/experimental/point_in_polygon.cuh>
-#include <cuspatial/experimental/type_utils.hpp>
 #include <cuspatial/vec_2d.hpp>
 
 #include <cudf/column/column.hpp>