diff --git a/CMakeLists.txt b/CMakeLists.txt index 33687cf..f547ef9 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -81,6 +81,7 @@ set( MAP_TYPE_SOURCES src/map_types/data_source_mapping.cpp src/map_types/dim_mapping.cpp src/map_types/expr_mapping.cpp + src/map_types/interp_mapping.cpp src/map_types/value_mapping.cpp ) @@ -116,6 +117,7 @@ set( MAP_TYPE_HEADERS src/map_types/data_source_mapping.hpp src/map_types/dim_mapping.hpp src/map_types/expr_mapping.hpp + src/map_types/interp_mapping.hpp src/map_types/map_arguments.hpp src/map_types/value_mapping.hpp ) diff --git a/include/libtokamap.hpp b/include/libtokamap.hpp index d382393..712d219 100644 --- a/include/libtokamap.hpp +++ b/include/libtokamap.hpp @@ -6,6 +6,7 @@ #include // IWYU pragma: export. #include // IWYU pragma: export. #include // IWYU pragma: export. +#include // IWYU pragma: export. #include // IWYU pragma: export. #include // IWYU pragma: export. #include // IWYU pragma: export. diff --git a/src/handlers/mapping_handler.cpp b/src/handlers/mapping_handler.cpp index d9a1ad7..0e0003d 100644 --- a/src/handlers/mapping_handler.cpp +++ b/src/handlers/mapping_handler.cpp @@ -30,6 +30,7 @@ #include "map_types/data_source_mapping.hpp" #include "map_types/dim_mapping.hpp" #include "map_types/expr_mapping.hpp" +#include "map_types/interp_mapping.hpp" #include "map_types/map_arguments.hpp" #include "map_types/value_mapping.hpp" #include "utils/algorithm.hpp" @@ -418,6 +419,31 @@ void init_expr_mapping(libtokamap::MappingStore& map_store, const libtokamap::Ma } } +void init_interp_mapping(libtokamap::MappingStore& map_store, const libtokamap::MappingName& mapping_name, + const nlohmann::json& value, libtokamap::MappingCounts& mapping_counts) +{ + for (const auto& required : {"input", "base", "target", "type"}) { + if (!value.contains(required)) { + throw libtokamap::ConfigurationError{"Required " + std::string{required} + + " argument not provided in INTERP mapping '" + mapping_name + "'"}; + } + } + + const auto input = value["input"].get(); + const auto base = value["base"].get(); + const auto target = value["target"].get(); + const auto interp_type = value["type"].get(); + + if (interp_type == libtokamap::InterpType::UNKNOWN) { + throw libtokamap::ConfigurationError{"Unknown interpolation type in INTERP mapping '" + mapping_name + "'"}; + } + + map_store.emplace(mapping_name, std::make_unique(input, base, target, interp_type)); + mapping_counts.increment(input); + mapping_counts.increment(base); + mapping_counts.increment(target); +} + void init_custom_mapping(libtokamap::MappingStore& map_store, const libtokamap::MappingName& mapping_name, const nlohmann::json& value, const std::vector& library_functions, libtokamap::MappingCounts& mapping_counts) @@ -480,6 +506,9 @@ libtokamap::MappingStore libtokamap::MappingHandler::init_mappings(const nlohman case MappingType::EXPR: init_expr_mapping(map_store, mapping_name, value, m_mapping_counts); break; + case MappingType::INTERP: + init_interp_mapping(map_store, mapping_name, value, m_mapping_counts); + break; case MappingType::CUSTOM: init_custom_mapping(map_store, mapping_name, value, m_library_functions, m_mapping_counts); break; diff --git a/src/map_types/base_mapping.hpp b/src/map_types/base_mapping.hpp index 15036e7..b3d35b0 100644 --- a/src/map_types/base_mapping.hpp +++ b/src/map_types/base_mapping.hpp @@ -9,13 +9,13 @@ namespace libtokamap { -enum class MappingType : uint8_t { UNKNOWN, VALUE, DATA_SOURCE, SLICE, EXPR, CUSTOM, DIM }; +enum class MappingType : uint8_t { UNKNOWN, VALUE, DATA_SOURCE, EXPR, CUSTOM, DIM, INTERP }; NLOHMANN_JSON_SERIALIZE_ENUM(MappingType, {{MappingType::UNKNOWN, ""}, // will default to this on no match {MappingType::VALUE, "VALUE"}, {MappingType::DATA_SOURCE, "DATA_SOURCE"}, - {MappingType::SLICE, "SLICE"}, {MappingType::EXPR, "EXPR"}, + {MappingType::INTERP, "INTERP"}, {MappingType::CUSTOM, "CUSTOM"}, {MappingType::DIM, "DIMENSION"}}) diff --git a/src/map_types/interp_mapping.cpp b/src/map_types/interp_mapping.cpp new file mode 100644 index 0000000..0da7f0c --- /dev/null +++ b/src/map_types/interp_mapping.cpp @@ -0,0 +1,184 @@ +#include "map_types/interp_mapping.hpp" + +#include +#include +#include +#include +#include +#include +#include +#include + +#include "exceptions/exceptions.hpp" +#include "map_types/map_arguments.hpp" +#include "utils/typed_data_array.hpp" + +namespace +{ + +using libtokamap::DataType; +using libtokamap::InterpType; +using libtokamap::TypedDataArray; + +std::vector to_doubles(const TypedDataArray& array, const std::string& name) +{ + if (array.data_type() == DataType::Float) { + const auto* data = array.data(); + return std::vector{data, data + array.size()}; + } + if (array.data_type() == DataType::Double) { + return array.to_vector(); + } + throw libtokamap::DataTypeError{"INTERP mapping '" + name + "' must be floating point, got " + + libtokamap::data_type_name(array.data_type())}; +} + +void check_increasing(const std::vector& base, const std::string& name) +{ + for (size_t idx = 1; idx < base.size(); ++idx) { + if (base[idx] <= base[idx - 1]) { + throw libtokamap::ProcessingError{"INTERP base '" + name + + "' must be strictly increasing, but is not at index " + + std::to_string(idx)}; + } + } +} + +/** + * @brief Indices of the two base points either side of the given point + * + * Points beyond either end of the base are clamped rather than extrapolated, + * and come back as a pair of equal indices for the caller to take as-is + */ +std::pair bracket(double point, const std::vector& base) +{ + // clamp, no extrapolate + if (point <= base.front()) { + return {0, 0}; + } + if (point >= base.back()) { + const size_t last = base.size() - 1; + return {last, last}; + } + + const auto upper = std::ranges::upper_bound(base, point); + const size_t high = static_cast(upper - base.begin()); + return {high - 1, high}; +} + +/** + * @brief Straight line between each pair of bracketing base points + * + * Every target point is an independent lookup into the base, so the targets + * need not be ordered + */ +std::vector interpolate_linear(const std::vector& target, const std::vector& base, + const std::vector& input) +{ + std::vector interpolated(target.size()); + std::ranges::transform(target, interpolated.begin(), [&](double point) { + const auto [low, high] = bracket(point, base); + if (low == high) { + return input[low]; + } + // fraction of the way between the bracketing base points + const double fraction = (point - base[low]) / (base[high] - base[low]); + return std::lerp(input[low], input[high], fraction); + }); + return interpolated; +} + +/** + * @brief Interpolate the input onto the target points, the base being the + * reference the input is defined against + * + * Each interpolation type owns its own loop, so that types needing a setup pass + * over the whole base (splines, for instance) have somewhere to do it once + * + * @param interp_type which interpolation to run, eg. linear + * @param target points to evaluate at, in any order + * @param base axis the input is sampled on, strictly increasing and the same + * length as the input + * @param input data values being interpolated, one per base point + * @return one interpolated value per target point, in target order + */ +std::vector interpolate(InterpType interp_type, const std::vector& target, + const std::vector& base, const std::vector& input) +{ + switch (interp_type) { + case InterpType::LINEAR: + return interpolate_linear(target, base, input); + case InterpType::UNKNOWN: + throw libtokamap::ProcessingError{"Unknown interpolation type"}; + // AJP: add more later + } + LIBTOKAMAP_UNREACHABLE +} + +} // namespace + +libtokamap::TypedDataArray libtokamap::InterpMapping::map_interp_args(const MapArguments& arguments, + const std::string& name) const +{ + if (!arguments.entries.contains(name)) { + throw libtokamap::MappingError{"Mapping '" + name + "' referenced by INTERP mapping not found"}; + } + return arguments.entries.at(name)->map(arguments); +} + +libtokamap::TypedDataArray libtokamap::InterpMapping::map(const MapArguments& arguments) const +{ + + const auto input_array = map_interp_args(arguments, m_input); + const auto base_array = map_interp_args(arguments, m_base); + const auto target_array = map_interp_args(arguments, m_target); + + // empty check + if (input_array.empty() || base_array.empty() || target_array.empty()) { + throw libtokamap::ProcessingError{"One (or more) of the INTERP parameters are empty"}; + } + + // rank 1D check + if (input_array.rank() != 1 || base_array.rank() != 1 || target_array.rank() != 1) { + throw libtokamap::ProcessingError{"Only 1D interpolation is supported" + " - please ensure all INTERP parameters are rank 1"}; + } + + // input size match == base size check + if (input_array.size() != base_array.size()) { + throw libtokamap::ProcessingError{"INTERP data '" + m_input + "' has " + std::to_string(input_array.size()) + + " elements but base '" + m_base + "' has " + + std::to_string(base_array.size())}; + } + + // Error for single point array check + if (base_array.size() < 2) { + throw libtokamap::ProcessingError{"INTERP base '" + m_base + + "' must have at least 2 points, got " + + std::to_string(base_array.size())}; + } + + const auto input = to_doubles(input_array, m_input); + const auto base = to_doubles(base_array, m_base); + const auto target = to_doubles(target_array, m_target); + + // ascending check + check_increasing(base, m_base); + + const auto interpolated = interpolate(m_interp_type, target, base, input); + + // use float if originally float + if (input_array.data_type() == DataType::Float) { + std::vector result; + result.reserve(interpolated.size()); + + std::transform( + interpolated.begin(), + interpolated.end(), + std::back_inserter(result), + [](double value) { return static_cast(value); }); + + return TypedDataArray{std::move(result)}; + } + return TypedDataArray{std::move(interpolated)}; +} diff --git a/src/map_types/interp_mapping.hpp b/src/map_types/interp_mapping.hpp new file mode 100644 index 0000000..1be029e --- /dev/null +++ b/src/map_types/interp_mapping.hpp @@ -0,0 +1,35 @@ +#pragma once + +#include + +#include "map_types/base_mapping.hpp" +#include "map_types/map_arguments.hpp" +#include "utils/typed_data_array.hpp" + +namespace libtokamap +{ + +enum class InterpType : short { UNKNOWN, LINEAR /*CUBIC, SPLINE, NEAREST*/ }; + +NLOHMANN_JSON_SERIALIZE_ENUM(InterpType, {{InterpType::UNKNOWN, ""}, // will default to this on no match + {InterpType::LINEAR, "LINEAR"}}) + +class InterpMapping : public Mapping +{ + public: + InterpMapping() = delete; + InterpMapping(std::string input, std::string base, std::string target, InterpType interp_type) + : m_input{std::move(input)}, m_base{std::move(base)}, m_target{std::move(target)}, m_interp_type{interp_type} {}; + + [[nodiscard]] TypedDataArray map(const MapArguments& arguments) const override; + + private: + std::string m_input; + std::string m_base; + std::string m_target; + InterpType m_interp_type; + + [[nodiscard]] TypedDataArray map_interp_args(const MapArguments& arguments, const std::string& name) const; +}; + +} // namespace libtokamap diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 618a7f6..ed1957d 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -14,6 +14,7 @@ set( TEST_SOURCES src/data_source_mapping_test.cpp src/dim_mapping_test.cpp src/indices_test.cpp + src/interp_mapping_test.cpp src/mapping_locator_test.cpp src/parse_slices_test.cpp src/render_test.cpp diff --git a/test/src/interp_mapping_test.cpp b/test/src/interp_mapping_test.cpp new file mode 100644 index 0000000..e47be61 --- /dev/null +++ b/test/src/interp_mapping_test.cpp @@ -0,0 +1,117 @@ +#include +#include +#include +#include +#include +#include + +#include "exceptions/exceptions.hpp" +#include "map_types/base_mapping.hpp" +#include "map_types/interp_mapping.hpp" +#include "map_types/map_arguments.hpp" +#include "map_types/value_mapping.hpp" +#include "utils/typed_data_array.hpp" + +#include +#include +#include + +using namespace libtokamap; + +namespace +{ + +using MapEntries = std::unordered_map>; + +MapArguments make_map_arguments(MapEntries& entries) +{ + static nlohmann::json empty_global_data = nlohmann::json::object(); + + constexpr bool trace_enabled = false; + constexpr bool cache_enabled = false; + constexpr RamCache* ram_cache = nullptr; + + return MapArguments(entries, empty_global_data, DataType::Float, 1, trace_enabled, cache_enabled, ram_cache); +} + +MapEntries make_entries(const nlohmann::json& signal, const nlohmann::json& time, const nlohmann::json& new_time) +{ + MapEntries entries; + entries["signal"] = std::make_unique(signal); + entries["time"] = std::make_unique(time); + entries["new_time"] = std::make_unique(new_time); + return entries; +} + +InterpMapping make_interp_mapping() +{ + return InterpMapping{"signal", "time", "new_time", InterpType::LINEAR}; +} + +} // namespace + +TEST_CASE("InterpMapping interpolates onto a new time base", "[interp_mapping]") +{ + // throw a few situations, on the nose, out of bounds, inbetween, all the rest + auto entries = make_entries( + {10.0, 20.0, 30.0, 40.0}, + {0.0, 1.0, 2.0, 3.0}, + {-1.0, 0.0, 0.5, 1.5, 2.0, 3.0, 5.0} + ); + MapArguments map_args = make_map_arguments(entries); + + const auto result = make_interp_mapping().map(map_args); + + REQUIRE(result.rank() == 1); + REQUIRE(result.size() == 7); + REQUIRE(result.data_type() == DataType::Float); + + // points outside the base are clamped rather than extrapolated + const std::vector expected{10.0F, 10.0F, 15.0F, 25.0F, 30.0F, 40.0F, 40.0F}; + const auto actual = result.to_vector(); + for (size_t idx = 0; idx < expected.size(); ++idx) { + INFO("index " << idx); + REQUIRE(actual[idx] == Catch::Approx(expected[idx])); + } +} + +TEST_CASE("InterpMapping rejects arguments it cannot interpolate", "[interp_mapping_errors]") +{ + SECTION("non floating point data") + { + auto entries = make_entries({10, 20}, {0.0, 1.0}, {0.5}); + MapArguments map_args = make_map_arguments(entries); + + REQUIRE_THROWS_AS(make_interp_mapping().map(map_args), DataTypeError); + REQUIRE_THROWS_WITH(make_interp_mapping().map(map_args), + Catch::Matchers::ContainsSubstring("must be floating point")); + } + + SECTION("input and base of different lengths") + { + auto entries = make_entries({10.0, 20.0, 30.0}, {0.0, 1.0}, {0.5}); + MapArguments map_args = make_map_arguments(entries); + + REQUIRE_THROWS_AS(make_interp_mapping().map(map_args), ProcessingError); + } + + SECTION("base that is not strictly increasing") + { + auto entries = make_entries({10.0, 20.0, 30.0, 40.0}, {0.0, 2.0, 1.0, 3.0}, {0.5}); + MapArguments map_args = make_map_arguments(entries); + + REQUIRE_THROWS_AS(make_interp_mapping().map(map_args), ProcessingError); + REQUIRE_THROWS_WITH(make_interp_mapping().map(map_args), + Catch::Matchers::ContainsSubstring("must be strictly increasing")); + } + + SECTION("a base with only one point") + { + auto entries = make_entries({10.0}, {0.0}, {0.5}); + MapArguments map_args = make_map_arguments(entries); + + REQUIRE_THROWS_AS(make_interp_mapping().map(map_args), ProcessingError); + REQUIRE_THROWS_WITH(make_interp_mapping().map(map_args), + Catch::Matchers::ContainsSubstring("must have at least 2 points")); + } +}