more helpers for NA bitmasking

kevinushey · kevinushey · commit fe6b87bbae15 · 2016-01-29T13:20:00.000-08:00
diff --git a/inst/examples/example-na-handling-variance.cpp b/inst/examples/example-na-handling-variance.cpp
@@ -9,28 +9,11 @@ using namespace RcppNT2;
 #include <Rcpp.h>
 using namespace Rcpp;
 
-struct IsNaN
-{
-  template <typename T>
-  T operator()(const T& data)
-  {
-    return nt2::if_else(data == data, bitwise::ones(), bitwise::zeroes());
-  }
-};
-
-struct NumNonZero
-{
-  template <typename T>
-  void operator()(const T& t)
-  {
-    result_ += nt2::sum(nt2::if_else(t != 0.0, 1.0, 0.0));
-  }
-
-  operator double() const { return result_; }
-
-  double result_ = 0.0;
-};
-
+// A 'mask-aware' Sum class. The call operator accepts both
+// 'data' and a 'mask' of the same time; one can use
+// 'nt2::bitwise_and' to apply the mask. In our case, we'll
+// have the mask set as bitwise '0' for NA', and bitwise '1'
+// otherwise.
 class Sum
 {
 public:
@@ -50,6 +33,9 @@ class Sum
 
 };
 
+// A 'mask-aware' Sum of Squares class. As above, the
+// call operator accepts both a 'data' and a 'mask'.
+// We use 'nt2::bitwise_and' to apply the mask.
 class SumOfSquares
 {
 public:
@@ -73,16 +59,14 @@ class SumOfSquares
 // [[Rcpp::export]]
 double simdVariance(NumericVector x)
 {
-  // Compute our NA bitmask. It's a vector o 'double's,
-  // as we want to ensure SIMD 'double' x 'double'
-  // instructions can be emitted when interacting with
-  // our bitmask. Note that we use our 'begin' and 'end'
-  // helpers
+  // Use helpers to compute the NA bitmask, as well as
+  // the number of non-NA elements in the vector.
   std::vector<double> naMask(x.size());
-  simdTransform(&x[0], &x[0] + x.size(), &naMask[0], IsNaN());
+  std::size_t n;
+  na::mask(pbegin(x), pend(x), pbegin(naMask), &n);
 
-  // Using our bitmask, compute the number of non-NA elements.
-  double N = simdFor(&naMask[0], &naMask[0] + naMask.size(), NumNonZero());
+  // Compute the number of non-NA elements.
+  std::size_t N = x.size() - n;
 
   // Compute the sum of our 'x' vector, discarding NAs.
   double total = simdFor(&x[0], &x[0] + x.size(), &naMask[0], Sum());
diff --git a/inst/include/RcppNT2.h b/inst/include/RcppNT2.h
@@ -44,6 +44,7 @@
 #include <RcppNT2/bitwise/bitwise.h>
 #include <RcppNT2/convert/convert.h>
 #include <RcppNT2/functor/functor.h>
+#include <RcppNT2/na/na.h>
 #include <RcppNT2/traits/traits.h>
 #include <RcppNT2/variadic/variadic.h>
 
diff --git a/inst/include/RcppNT2/core/core.h b/inst/include/RcppNT2/core/core.h
@@ -6,13 +6,13 @@
 namespace RcppNT2 {
 
 template <typename T>
-auto begin(T&& t) -> decltype(&std::forward<T>(t)[0])
+auto pbegin(T&& t) -> decltype(&std::forward<T>(t)[0])
 {
   return &std::forward<T>(t)[0];
 }
 
 template <typename T>
-auto end(T&& t) -> decltype(&std::forward<T>(t)[0])
+auto pend(T&& t) -> decltype(&std::forward<T>(t)[0])
 {
   return &std::forward<T>(t)[0] + std::forward<T>(t).size();
 }
diff --git a/inst/include/RcppNT2/na/na.h b/inst/include/RcppNT2/na/na.h
@@ -0,0 +1,78 @@
+#ifndef RCPP_NT2_NA_NA_H
+#define RCPP_NT2_NA_NA_H
+
+#include <RcppNT2/bitwise/bitwise.h>
+
+namespace RcppNT2 {
+namespace na {
+
+namespace detail {
+
+struct Masker
+{
+  template <typename T>
+  T operator()(const T& data)
+  {
+    static const double ones   = bitwise::ones();
+    static const double zeroes = bitwise::zeroes();
+    return nt2::if_else(data == data, ones, zeroes);
+  }
+};
+
+struct MaskerWithCount
+{
+
+  template <typename T>
+  T operator()(const T& data)
+  {
+    static const double ones   = bitwise::ones();
+    static const double zeroes = bitwise::zeroes();
+    auto&& compare = data == data;
+    count_ += nt2::sum(nt2::if_else(compare, 1.0, 0.0));
+    return nt2::if_else(compare, ones, zeroes);
+  }
+
+  std::size_t count() { return count_; }
+  std::size_t count_ = 0;
+};
+
+struct Counter
+{
+  template <typename T>
+  void operator()(const T& t)
+  {
+    count_ += nt2::sum(t != t);
+  }
+
+  template <typename T>
+  operator T() const { return T(count_); }
+
+  std::size_t count_ = 0;
+};
+
+} // namespace detail
+
+template <typename T, typename U>
+void mask(const T* begin, const T* end, U* out)
+{
+  simdTransform(begin, end, out, detail::Masker());
+}
+
+template <typename T, typename U, typename V>
+void mask(const T* begin, const T* end, U* out, V* count)
+{
+  detail::MaskerWithCount masker;
+  simdTransform(begin, end, out, std::reference_wrapper<detail::MaskerWithCount>(masker));
+  *count = masker.count();
+}
+
+template <typename T>
+std::size_t count(const T& t)
+{
+  return simdReduce(pbegin(t), pend(t), std::size_t(0), detail::Counter());
+}
+
+} // namespace na
+} // namespace RcppNT2
+
+#endif /* RCPP_NT2_NA_NA_H */

Original file line number	Diff line number	Diff line change
`@@ -6,13 +6,13 @@`
`6`	`6`	`namespace RcppNT2 {`
`7`	`7`
`8`	`8`	`template <typename T>`
`9`		`-auto begin(T&& t) -> decltype(&std::forward<T>(t)[0])`
	`9`	`+auto pbegin(T&& t) -> decltype(&std::forward<T>(t)[0])`
`10`	`10`	`{`
`11`	`11`	`return &std::forward<T>(t)[0];`
`12`	`12`	`}`
`13`	`13`
`14`	`14`	`template <typename T>`
`15`		`-auto end(T&& t) -> decltype(&std::forward<T>(t)[0])`
	`15`	`+auto pend(T&& t) -> decltype(&std::forward<T>(t)[0])`
`16`	`16`	`{`
`17`	`17`	`return &std::forward<T>(t)[0] + std::forward<T>(t).size();`
`18`	`18`	`}`