From ea5bbc18953ae585377e975ff8700d2486e3bd66 Mon Sep 17 00:00:00 2001 From: Daniel Baston Date: Tue, 19 Dec 2023 20:05:57 -0500 Subject: [PATCH] RasterSource: Add read_empty for use in dispatch --- src/operation.cpp | 89 ++++++++++++++++++++++++++++++++++++++++++ src/operation.h | 88 +---------------------------------------- src/raster_source.h | 14 ++++++- src/stats_registry.cpp | 8 ++-- 4 files changed, 106 insertions(+), 93 deletions(-) diff --git a/src/operation.cpp b/src/operation.cpp index b4ff0f6c..51ac0de6 100644 --- a/src/operation.cpp +++ b/src/operation.cpp @@ -4,6 +4,95 @@ namespace exactextract { +void +Operation::set_result(const StatsRegistry& reg, const Feature& f_in, Feature& f_out) const +{ + static const StatsRegistry::RasterStatsVariant empty_stats = RasterStats(); + + constexpr bool write_if_missing = true; // should we set attribute values if the feature did not intersect the raster? + if (!write_if_missing && !reg.contains(f_in, *this)) { + return; + } + + const auto& stats = reg.contains(f_in, *this) ? reg.stats(f_in, *this) : empty_stats; + + // Construct + const auto& empty_rast = values->read_empty(); + + using missing_value_t = std::variant; + auto missing = std::visit([](const auto& r) -> missing_value_t { + if (r->has_nodata()) { + return r->nodata(); + } + return std::numeric_limits::quiet_NaN(); + }, + empty_rast); + + if (stat == "mean") { + std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.mean()); }, stats); + } else if (stat == "sum") { + std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.sum()); }, stats); + } else if (stat == "count") { + std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.count()); }, stats); + } else if (stat == "weighted_mean") { + std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.weighted_mean()); }, stats); + } else if (stat == "weighted_sum") { + std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.weighted_sum()); }, stats); + } else if (stat == "min") { + std::visit([&f_out, this](const auto& x, const auto& m) { f_out.set(m_field_names[0], x.min().value_or(m)); }, stats, missing); + } else if (stat == "max") { + std::visit([&f_out, this](const auto& x, const auto& m) { f_out.set(m_field_names[0], x.max().value_or(m)); }, stats, missing); + } else if (stat == "majority" || stat == "mode") { + std::visit([&f_out, this](const auto& x, const auto& m) { f_out.set(m_field_names[0], x.mode().value_or(m)); }, stats, missing); + } else if (stat == "minority") { + std::visit([&f_out, this](const auto& x, const auto& m) { f_out.set(m_field_names[0], x.minority().value_or(m)); }, stats, missing); + } else if (stat == "variety") { + std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.variety()); }, stats); + } else if (stat == "stdev") { + std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.stdev()); }, stats); + } else if (stat == "weighted_stdev") { + std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.weighted_stdev()); }, stats); + } else if (stat == "variance") { + std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.variance()); }, stats); + } else if (stat == "weighted_variance") { + std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.weighted_variance()); }, stats); + } else if (stat == "coefficient_of_variation") { + std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.coefficient_of_variation()); }, stats); + } else if (stat == "median") { + std::visit([&f_out, this](const auto& x, const auto& m) { f_out.set(m_field_names[0], x.quantile(0.5).value_or(m)); }, stats, missing); + } else if (stat == "quantile") { + std::visit([&f_out, this](const auto& x, const auto& m) { + for (std::size_t i = 0; i < m_quantiles.size(); i++) { + f_out.set(m_field_names[i], x.quantile(m_quantiles[i]).value_or(m)); + } + }, + stats, + missing); + } else if (stat == "frac") { + std::visit([&f_out](const auto& x) { + for (const auto& value : x) { + std::stringstream s; + s << "frac_" << value; + + f_out.set(s.str(), x.frac(value).value_or(0)); + } + }, + stats); + } else if (stat == "weighted_frac") { + std::visit([&f_out](const auto& x) { + for (const auto& value : x) { + std::stringstream s; + s << "weighted_frac_" << value; + + f_out.set(s.str(), x.weighted_frac(value).value_or(0)); + } + }, + stats); + } else { + throw std::runtime_error("Unhandled stat: " + stat); + } +} + void Operation::setQuantileFieldNames() { diff --git a/src/operation.h b/src/operation.h index 8a0b6b96..41d8b1c8 100644 --- a/src/operation.h +++ b/src/operation.h @@ -107,93 +107,7 @@ class Operation setQuantileFieldNames(); } - virtual void set_result(const StatsRegistry& reg, const Feature& f_in, Feature& f_out) const - { - static const StatsRegistry::RasterStatsVariant empty_stats = RasterStats(); - - constexpr bool write_if_missing = true; // should we set attribute values if the feature did not intersect the raster? - if (!write_if_missing && !reg.contains(f_in, *this)) { - return; - } - - const auto& stats = reg.contains(f_in, *this) ? reg.stats(f_in, *this) : empty_stats; - - // FIXME don't read an empty box every time, maybe cache this in the source? - auto empty_rast = values->read_box(Box::make_empty()); - - auto missing = std::visit([](const auto& r) { - std::variant ret = std::numeric_limits::quiet_NaN(); - if (r->has_nodata()) { - ret = r->nodata(); - } - return ret; - }, - empty_rast); - - if (stat == "mean") { - std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.mean()); }, stats); - } else if (stat == "sum") { - std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.sum()); }, stats); - } else if (stat == "count") { - std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.count()); }, stats); - } else if (stat == "weighted_mean") { - std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.weighted_mean()); }, stats); - } else if (stat == "weighted_sum") { - std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.weighted_sum()); }, stats); - } else if (stat == "min") { - std::visit([&f_out, this](const auto& x, const auto& m) { f_out.set(m_field_names[0], x.min().value_or(m)); }, stats, missing); - } else if (stat == "max") { - std::visit([&f_out, this](const auto& x, const auto& m) { f_out.set(m_field_names[0], x.max().value_or(m)); }, stats, missing); - } else if (stat == "majority" || stat == "mode") { - std::visit([&f_out, this](const auto& x, const auto& m) { f_out.set(m_field_names[0], x.mode().value_or(m)); }, stats, missing); - } else if (stat == "minority") { - std::visit([&f_out, this](const auto& x, const auto& m) { f_out.set(m_field_names[0], x.minority().value_or(m)); }, stats, missing); - } else if (stat == "variety") { - std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.variety()); }, stats); - } else if (stat == "stdev") { - std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.stdev()); }, stats); - } else if (stat == "weighted_stdev") { - std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.weighted_stdev()); }, stats); - } else if (stat == "variance") { - std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.variance()); }, stats); - } else if (stat == "weighted_variance") { - std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.weighted_variance()); }, stats); - } else if (stat == "coefficient_of_variation") { - std::visit([&f_out, this](const auto& x) { f_out.set(m_field_names[0], x.coefficient_of_variation()); }, stats); - } else if (stat == "median") { - std::visit([&f_out, this](const auto& x, const auto& m) { f_out.set(m_field_names[0], x.quantile(0.5).value_or(m)); }, stats, missing); - } else if (stat == "quantile") { - std::visit([&f_out, this](const auto& x, const auto& m) { - for (std::size_t i = 0; i < m_quantiles.size(); i++) { - f_out.set(m_field_names[i], x.quantile(m_quantiles[i]).value_or(m)); - } - }, - stats, - missing); - } else if (stat == "frac") { - std::visit([&f_out](const auto& x) { - for (const auto& value : x) { - std::stringstream s; - s << "frac_" << value; - - f_out.set(s.str(), x.frac(value).value_or(0)); - } - }, - stats); - } else if (stat == "weighted_frac") { - std::visit([&f_out](const auto& x) { - for (const auto& value : x) { - std::stringstream s; - s << "weighted_frac_" << value; - - f_out.set(s.str(), x.weighted_frac(value).value_or(0)); - } - }, - stats); - } else { - throw std::runtime_error("Unhandled stat: " + stat); - } - } + virtual void set_result(const StatsRegistry& reg, const Feature& f_in, Feature& f_out) const; std::string stat; std::string name; diff --git a/src/raster_source.h b/src/raster_source.h index a7c7c7f4..a459a134 100644 --- a/src/raster_source.h +++ b/src/raster_source.h @@ -1,4 +1,4 @@ -// Copyright (c) 2020 ISciences, LLC. +// Copyright (c) 2020-2023 ISciences, LLC. // All rights reserved. // // This software is licensed under the Apache License, Version 2.0 (the "License"). @@ -27,6 +27,15 @@ class RasterSource virtual const Grid& grid() const = 0; virtual RasterVariant read_box(const Box& box) = 0; + const RasterVariant& read_empty() + { + if (!m_empty) { + m_empty = std::make_unique(read_box(Box::make_empty())); + } + + return *m_empty; + } + virtual ~RasterSource() = default; void set_name(const std::string& name) @@ -34,13 +43,14 @@ class RasterSource m_name = name; } - std::string name() const + const std::string& name() const { return m_name; } private: std::string m_name; + std::unique_ptr m_empty; }; } diff --git a/src/stats_registry.cpp b/src/stats_registry.cpp index e81dc2f8..0c3d903e 100644 --- a/src/stats_registry.cpp +++ b/src/stats_registry.cpp @@ -52,13 +52,13 @@ StatsRegistry::stats(const Feature& feature, const Operation& op, bool store_val auto it = stats_for_feature.find(op.key()); if (it == stats_for_feature.end()) { - // FIXME come up with a more direct way to probe the type of the raster - auto rast = op.values->read_box(Box::make_empty()); + // Construct a RasterStats with the correct value type for this Operation. + const auto& rast = op.values->read_empty(); it = std::visit([&stats_for_feature, &op, store_values](const auto& r) { - using rtype = typename std::remove_reference_t::value_type; + using value_type = typename std::remove_reference_t::value_type; - return stats_for_feature.emplace(op.key(), RasterStats(store_values)); + return stats_for_feature.emplace(op.key(), RasterStats(store_values)); }, rast) .first;