From 3d27e6e27ac119ff734b38f30569db8665156685 Mon Sep 17 00:00:00 2001 From: AlePalu Date: Mon, 25 Nov 2024 20:51:51 +0100 Subject: [PATCH] geoframe enables indexed access to layers, as plain_data_layer instances. Iterator on layers --- fdaPDE/geoframe/areal_layer.h | 1 + fdaPDE/geoframe/data_layer.h | 1 + fdaPDE/geoframe/geoframe.h | 111 +++++++++++++++++++++++++++++----- fdaPDE/geoframe/point_layer.h | 2 +- 4 files changed, 100 insertions(+), 15 deletions(-) diff --git a/fdaPDE/geoframe/areal_layer.h b/fdaPDE/geoframe/areal_layer.h index f0aaac4..e50336c 100644 --- a/fdaPDE/geoframe/areal_layer.h +++ b/fdaPDE/geoframe/areal_layer.h @@ -63,6 +63,7 @@ template struct areal_layer { data_ = storage_t(std::forward(data)...); } storage_t& data() { return data_; } + const storage_t& data() const { return data_; } // geometry const MultiPolygon& geometry(int i) const { return regions_->operator[](i); } // computes measures of subdomains diff --git a/fdaPDE/geoframe/data_layer.h b/fdaPDE/geoframe/data_layer.h index aa6de2f..2a8d865 100644 --- a/fdaPDE/geoframe/data_layer.h +++ b/fdaPDE/geoframe/data_layer.h @@ -179,6 +179,7 @@ template struct plain_col_view { size_t cols() const { return 1; } size_t size() const { return data_->rows(); } index_t id() const { return col_; } + const storage_t& data() const { return slice_; } const auto& field_descriptor() const { return data_->field_descriptor(colname_); } const std::string& colname() const { return colname_; } internals::dtype type_id() const { return type_id_; } diff --git a/fdaPDE/geoframe/geoframe.h b/fdaPDE/geoframe/geoframe.h index 8aca208..fef8b72 100644 --- a/fdaPDE/geoframe/geoframe.h +++ b/fdaPDE/geoframe/geoframe.h @@ -40,7 +40,9 @@ template struct GeoFrame { fdapde_static_assert(Order_ > 1, GEOFRAME_MUST_HAVE_ORDER_TWO_OR_HIGHER); private: using This = GeoFrame; - using layers_t = std::tuple, internals::areal_layer>; + using point_layer_t = internals::point_layer; + using areal_layer_t = internals::areal_layer; + using layers_t = std::tuple; template using LayerMap_ = std::tuple...>; template auto& fetch_(U& u) { return std::get::index>(u); } template const auto& fetch_(const U& u) const { @@ -64,16 +66,18 @@ template struct GeoFrame { using Triangulation = Triangulation_; // constructors - GeoFrame() noexcept : triangulation_(nullptr), layers_() { } + GeoFrame() noexcept : triangulation_(nullptr), layers_(), n_layers_(0) { } explicit GeoFrame(Triangulation_& triangulation) noexcept : - triangulation_(std::addressof(triangulation)), layers_() { } - + triangulation_(std::addressof(triangulation)), layers_(), n_layers_(0) { } + // modifiers // multipoint layer with geometrical locations at mesh nodes void push(const std::string& name, layer_t::point_t) { geoframe_assert(!name.empty() && !has_layer(name), "empty or duplicated name."); using layer_t = internals::point_layer; fetch_(layers_).insert({name, layer_t(name, this)}); + idx_to_layer_name_[n_layers_] = name; + n_layers_++; return; } template @@ -107,6 +111,8 @@ template struct GeoFrame { std::make_shared>(Eigen::Map>(coords.data(), n_rows, n_cols)); fetch_(layers_).insert({name, layer_t(name, this, coords_ptr)}); } + idx_to_layer_name_[n_layers_] = name; + n_layers_++; return; } // packed multipoint layers sharing the same locations @@ -134,6 +140,8 @@ template struct GeoFrame { std::string name(*it); geoframe_assert(!name.empty() && !has_layer(name), "empty or duplicated name."); fetch_(layers_).insert({name, layer_t(name, this, coords_ptr)}); + idx_to_layer_name_[n_layers_] = name; + n_layers_++; } } @@ -145,6 +153,8 @@ template struct GeoFrame { using areal_t = std::vector>; std::shared_ptr regions_ptr = std::make_shared(regions); fetch_(layers_).insert({name, layer_t(name, this, regions_ptr)}); + idx_to_layer_name_[n_layers_] = name; + n_layers_++; return; } // packed areal layers sharing the same regions @@ -158,8 +168,11 @@ template struct GeoFrame { // allocate shared memory std::shared_ptr regions_ptr = std::make_shared(regions); for (auto it = layers.begin(); it != layers.end(); ++it) { - geoframe_assert(!it->empty() && !has_layer(*it), "empty or duplicated name."); - fetch_(layers_).insert({*it, layer_t(*it, this, regions_ptr)}); + std::string name(*it); + geoframe_assert(!it->empty() && !has_layer(name), "empty or duplicated name."); + fetch_(layers_).insert({name, layer_t(name, this, regions_ptr)}); + idx_to_layer_name_[n_layers_] = name; + n_layers_++; } return; } @@ -208,6 +221,8 @@ template struct GeoFrame { void push(const std::string& name, Layer&& layer) { using layer_t = std::decay_t; fetch_(layers_).insert({name, layer}); + idx_to_layer_name_[n_layers_] = name; + n_layers_++; return; } void erase(const std::string& layer_name) { @@ -224,6 +239,17 @@ template struct GeoFrame { ...); }, layers_); + // update idx - layer_name mapping + int i = 0; + for (; i < n_layers_; ++i) { + if (idx_to_layer_name_.at(i) == layer_name) { + idx_to_layer_name_.erase(i); + break; + } + } + for (; i < n_layers_ - 1; ++i) { idx_to_layer_name_[i] = idx_to_layer_name_.at(i + 1); } + idx_to_layer_name_.erase(n_layers_ - 1); + n_layers_--; return; } @@ -238,6 +264,8 @@ template struct GeoFrame { // move data in memory buffer get_as(layer_t::point, name).data().assign_inplace_from(csv.data()); get_as(layer_t::point, name).set_colnames(csv.colnames()); + idx_to_layer_name_[n_layers_] = name; + n_layers_++; return; } template @@ -254,8 +282,11 @@ template struct GeoFrame { // move data in memory buffer get_as(layer_t::point, name).data().assign_inplace_from(csv.data()); get_as(layer_t::point, name).set_colnames(csv.colnames()); + idx_to_layer_name_[n_layers_] = name; + n_layers_++; return; } + void load_shp(const std::string& name, const std::string& file_name) { std::string file_name_ = std::filesystem::current_path().string() + "/" + file_name; geoframe_assert(std::filesystem::exists(file_name_), "file " + file_name_ + " not found."); @@ -274,25 +305,27 @@ template struct GeoFrame { std::vector>> int_data; std::vector>> dbl_data; std::vector>> str_data; - //TODO: std::vector>> bin_data; + // TODO: std::vector>> bin_data; for (const auto& [name, field_type] : shp.field_descriptors()) { if (field_type == 'N') { dbl_data.emplace_back(name, shp.get(name)); } if (field_type == 'C') { str_data.emplace_back(name, shp.get(name)); } } get_as(layer_t::areal, name).set_data(int_data, dbl_data, str_data/*, bin_data*/); - + idx_to_layer_name_[n_layers_] = name; + n_layers_++; + // TODO: we should reorder the fields in the same order they come from the shp } - // iterators - auto begin(layer_t::point_t) { return fetch_>(layers_).begin(); } - auto end(layer_t::point_t) { return fetch_>(layers_).end(); } // layer access template auto& get_as(Tag, const std::string& name) { return get_as_>(name); } + template auto& get_as(Tag t, int idx) { return get_as(t, idx_to_layer_name_.at(idx)); } template const auto& get_as(Tag, const std::string& name) const { return get_as_>(name); } + template const auto& get_as(Tag t, int idx) const { return get_as(t, idx_to_layer_name_.at(idx)); } + // observers bool has_layer(const std::string& name) const { // search for layer_name in each layer type bool found_ = false; @@ -312,12 +345,60 @@ template struct GeoFrame { layers_); return found_; } - std::optional layer_type(const std::string& name) const { + std::optional layer_category(const std::string& name) const { geoframe_assert(has_layer(name), std::string("key " + name + " not found.")); - if (fetch_>(layers_).contains(name)) return internals::ltype::point; - if (fetch_>(layers_).contains(name)) return internals::ltype::areal; + if (fetch_(layers_).contains(name)) return internals::ltype::point; + if (fetch_(layers_).contains(name)) return internals::ltype::areal; return std::nullopt; } + std::optional layer_category(int idx) const { return layer_category(idx_to_layer_name_.at(idx)); } + int n_layers() const { return n_layers_; } + // indexed access + const internals::plain_data_layer& operator[](int idx) const { + geoframe_assert(idx < n_layers_, "out of bound access."); + auto get_plain_data_ = + [&, this]([[maybe_unused]] LayerType l, const std::string& layer_name) -> const void* { + if (fetch_(layers_).contains(layer_name)) { + return std::addressof(fetch_(layers_).at(layer_name).data()); + } + return nullptr; + }; + std::string name = idx_to_layer_name_.at(idx); + const void* layer_ptr; + // as idx < n_layers_, it is guaranteed that one of the branches will be taken + if (fetch_(layers_).contains(name)) layer_ptr = get_plain_data_(point_layer_t {}, name); + if (fetch_(layers_).contains(name)) layer_ptr = get_plain_data_(areal_layer_t {}, name); + return *reinterpret_cast*>(layer_ptr); + } + // iterator + class iterator { + const GeoFrame* gf_; + int index_; + public: + using value_type = internals::plain_data_layer; + using pointer = std::add_pointer_t; + using reference = std::add_lvalue_reference_t; + using size_type = std::size_t; + using difference_type = std::ptrdiff_t; + using iterator_category = std::forward_iterator_tag; + + iterator(const GeoFrame* gf, int index) : gf_(gf), index_(index) { } + const value_type* operator->() const { return std::addressof(gf_->operator[](index_)); } + const value_type& operator*() const { return gf_->operator[](index_); } + iterator& operator++() { + index_++; + return *this; + } + internals::ltype category() const { return *(gf_->layer_category(index_)); } + const reference data() const { return operator*(); } + template const auto& as(Tag t) const { + return gf_->get_as(t, gf_->idx_to_layer_name_.at(index_)); + } + friend bool operator!=(const iterator& lhs, const iterator& rhs) { return lhs.index_ != rhs.index_; } + friend bool operator==(const iterator& lhs, const iterator& rhs) { return lhs.index_ == rhs.index_; } + }; + iterator begin() const { return iterator(this, 0); } + iterator end() const { return iterator(this, n_layers_); } // geometry access Triangulation_& triangulation() { return *triangulation_; } const Triangulation_& triangulation() const { return *triangulation_; } @@ -344,6 +425,8 @@ template struct GeoFrame { // data members Triangulation_* triangulation_ = nullptr; LayerMap layers_ {}; + int n_layers_ = 0; + std::unordered_map idx_to_layer_name_; }; // // regions is a vector with the same number of elements as number of cells diff --git a/fdaPDE/geoframe/point_layer.h b/fdaPDE/geoframe/point_layer.h index 1b35f46..6f96bca 100644 --- a/fdaPDE/geoframe/point_layer.h +++ b/fdaPDE/geoframe/point_layer.h @@ -61,7 +61,7 @@ template struct point_layer { } storage_t& data() { return data_; } const storage_t& data() const { return data_; } - + // geometry const DMatrix& coordinates() const { return coords_ ? *coords_ : triangulation().nodes(); } Triangulation& triangulation() { return geoframe_->triangulation(); } const Triangulation& triangulation() const { return geoframe_->triangulation(); }