compute_weights_gauss: call compute_weights_dist() if weights_sum<1e-9

src/Remapgrid.cc

@@ -344,7 +344,31 @@ get_parameter()
 }
 
 static void
-print_parameter(const RemapgridParams &params)
+verify_knn_parameter(KnnParams const &knnParams)
+{
+  if (knnParams.weighted == WeightingMethod::gaussWeighted)
+    {
+      if (knnParams.k <= 1) cdo_abort("Parameter k must be greater than 1!");
+      if (knnParams.kMin <= 1) cdo_abort("Parameter kmin must be greater than 1!");
+    }
+}
+
+static void
+print_knn_parameter(KnnParams const &knnParams)
+{
+  std::stringstream outbuffer;
+  outbuffer << "KNN parameter: ";
+  outbuffer << "k=" << knnParams.k;
+  outbuffer << ", kmin=" << knnParams.kMin;
+  outbuffer << ", weighted=" << weightingMethod_to_string(knnParams.weighted);
+  outbuffer << ", gauss_scale=" << knnParams.gaussScale;
+  outbuffer << ", extrapolate=" << knnParams.extrapolate;
+
+  cdo_verbose("%s", outbuffer.str());
+}
+
+static void
+print_parameter(RemapgridParams const &params)
 {
   std::stringstream outbuffer;
   outbuffer << "grid=" << params.gridString;
@@ -458,6 +482,8 @@ public:
             targetGridName = remapParams.gridString;
             knnParams = remapParams.knnParams;
             if (knnParams.kMin == 0) knnParams.kMin = knnParams.k;
+            verify_knn_parameter(knnParams);
+            print_knn_parameter(knnParams);
           }
         else
           {

src/knndata.cc

@@ -134,7 +134,8 @@ KnnData::compute_weights_gauss()
 {
   double gauss_scale = m_gaussScale;
   auto n = m_numNeighbors;
-  for (size_t i = 0; i < n; ++i) { m_dist[i] = m_dist[i] * m_dist[i]; }
+  std::vector<double> weights(n);
+  for (size_t i = 0; i < n; ++i) { weights[i] = m_dist[i] * m_dist[i]; }
 
   // a) compute sum of source point distances
   double src_distances_sum = 0.0;
@@ -152,8 +153,8 @@ KnnData::compute_weights_gauss()
   double weights_sum = 0.0;
   for (size_t i = 0; i < n; ++i)
     {
-      m_dist[i] = exp(scale * m_dist[i]);
-      weights_sum += m_dist[i];
+      weights[i] = exp(scale * weights[i]);
+      weights_sum += weights[i];
     }
 
   // If the sum of the weights is very low, which can happen in case
@@ -164,12 +165,12 @@ KnnData::compute_weights_gauss()
       // point cannot be computed. Therefore, the normalisation would
       // generate NaN's. Hence we fall back to inverse distance weighted
       // averge for this target point.
-      // compute_weights_dist(m_tgtCoord, m_srcCoords, n, m_dist, gauss_scale);  <<<<<---------------
+      compute_weights_dist();
       return n;
     }
 
   // d) scale weights such that SUM(w_i) == 1
-  for (size_t i = 0; i < n; ++i) m_dist[i] /= weights_sum;
+  for (size_t i = 0; i < n; ++i) m_dist[i] = weights[i] / weights_sum;
 
   return n;
 }