range.hpp

Go to the documentation of this file.

/* Copyright (c) 2025 Otto Link. Distributed under the terms of the GNU General
   Public License. The full license is in the file LICENSE, distributed with
   this software. */
#pragma once
#include <fstream>
#include <sstream>
#include <string>
#include <vector>
 
#include "point_sampler/point.hpp"
 
namespace ps
{
 
template <typename T, size_t N>
std::vector<Point<T, N>> filter_points_in_range(
    const std::vector<Point<T, N>>       &points,
    const std::array<std::pair<T, T>, N> &axis_ranges)
{
  std::vector<Point<T, N>> filtered;
  filtered.reserve(points.size());
 
  for (const auto &p : points)
  {
    bool inside = true;
    for (size_t i = 0; i < N; ++i)
    {
      const auto &[min_val, max_val] = axis_ranges[i];
      if (p[i] < min_val || p[i] > max_val)
      {
        inside = false;
        break;
      }
    }
    if (inside)
      filtered.push_back(p);
  }
 
  return filtered;
}
 
template <typename T, std::size_t N, typename Func>
std::vector<Point<T, N>> filter_points_function(const std::vector<Point<T, N>> &points,
                                                Func                            fn)
{
  std::vector<Point<T, N>> filtered;
  filtered.reserve(points.size());
 
  for (const auto &p : points)
  {
    if (fn(p) != T(0)) // keep point if function value is not zero
      filtered.push_back(p);
  }
 
  return filtered;
}
 
template <typename T, size_t N>
void refit_points_to_range(std::vector<Point<T, N>>             &points,
                           const std::array<std::pair<T, T>, N> &target_ranges)
{
  if (points.empty())
    return;
 
  std::array<T, N> min_vals, max_vals;
 
  // Initialize min/max
  for (size_t d = 0; d < N; ++d)
  {
    min_vals[d] = points[0][d];
    max_vals[d] = points[0][d];
  }
 
  // Compute bounding box
  for (const auto &p : points)
  {
    for (size_t d = 0; d < N; ++d)
    {
      min_vals[d] = std::min(min_vals[d], p[d]);
      max_vals[d] = std::max(max_vals[d], p[d]);
    }
  }
 
  // Apply linear mapping to each point
  for (auto &p : points)
  {
    for (size_t d = 0; d < N; ++d)
    {
      const T in_min = min_vals[d];
      const T in_max = max_vals[d];
      const T out_min = target_ranges[d].first;
      const T out_max = target_ranges[d].second;
 
      if (std::abs(in_max - in_min) < T(1e-12))
      {
        // Degenerate axis: center the target range
        p[d] = (out_min + out_max) / T(2);
      }
      else
      {
        const T t = (p[d] - in_min) / (in_max - in_min);
        p[d] = out_min + t * (out_max - out_min);
      }
    }
  }
}
 
template <typename T, size_t N>
void rescale_points(std::vector<Point<T, N>>             &points,
                    const std::array<std::pair<T, T>, N> &ranges)
{
  for (auto &pt : points)
    for (size_t d = 0; d < N; ++d)
      pt[d] = ranges[d].first + pt[d] * (ranges[d].second - ranges[d].first);
}
 
} // namespace ps