array_weighted_mean: added OpenMP pragma.

src/varray.cc

@@ -667,6 +667,11 @@ template <typename T>
 double
 varray_weighted_mean(const size_t len, const Varray<T> &v, const Varray<double> &w, const T missval)
 {
+  auto f_weighted_mean = [](const auto aw, const auto a, auto &sum, auto &sumw) {
+    sum += aw * a;
+    sumw += aw;
+  };
+
   assert(len > 0);
   assert(v.size() > 0);
   assert(len <= v.size());
@@ -674,13 +679,12 @@ varray_weighted_mean(const size_t len, const Varray<T> &v, const Varray<double>
 
   double sum = 0.0, sumw = 0.0;
 
-  for (size_t i = 0; i < len; ++i)
-    {
-      sum += w[i] * v[i];
-      sumw += w[i];
-    }
+#ifdef HAVE_OPENMP4
+#pragma omp parallel for default(shared) schedule(static) reduction(+ : sum,sumw)
+#endif
+  for (size_t i = 0; i < len; ++i) f_weighted_mean(w[i], v[i], sum, sumw);
 
-  return IS_EQUAL(sumw, 0.) ? missval : sum / sumw;
+  return is_equal(sumw, 0.0) ? missval : sum / sumw;
 }
 
 // Explicit instantiation
@@ -691,6 +695,14 @@ template <typename T>
 double
 varray_weighted_mean_mv(const size_t len, const Varray<T> &v, const Varray<double> &w, const T missval)
 {
+  auto f_weighted_mean_mv = [](const auto aw, const auto a, const auto mv_a, auto &sum, auto &sumw, auto is_EQ) {
+    if (!is_EQ(a, mv_a) && !is_EQ(aw, mv_a))
+      {
+        sum += aw * a;
+        sumw += aw;
+      }
+  };
+
   assert(len > 0);
   assert(v.size() > 0);
   assert(len <= v.size());
@@ -699,12 +711,20 @@ varray_weighted_mean_mv(const size_t len, const Varray<T> &v, const Varray<doubl
   const double missval1 = missval, missval2 = missval;
   double sum = 0.0, sumw = 0.0;
 
-  for (size_t i = 0; i < len; ++i)
-    if (!DBL_IS_EQUAL(v[i], missval1) && !DBL_IS_EQUAL(w[i], missval1))
-      {
-        sum += w[i] * v[i];
-        sumw += w[i];
-      }
+  if (std::isnan(missval))
+    {
+#ifdef HAVE_OPENMP4
+#pragma omp parallel for default(shared) schedule(static) reduction(+ : sum,sumw)
+#endif
+      for (size_t i = 0; i < len; ++i) f_weighted_mean_mv(w[i], v[i], missval1, sum, sumw, dbl_is_equal);
+    }
+  else
+    {
+#ifdef HAVE_OPENMP4
+#pragma omp parallel for default(shared) schedule(static) reduction(+ : sum,sumw)
+#endif
+      for (size_t i = 0; i < len; ++i) f_weighted_mean_mv(w[i], v[i], missval1, sum, sumw, is_equal);
+    }
 
   return DIVMN(sum, sumw);
 }