Refactor

src/libhictk/file/include/hictk/file.hpp

@@ -14,6 +14,7 @@
 #include <vector>
 
 #include "hictk/balancing/methods.hpp"
+#include "hictk/balancing/weights.hpp"
 #include "hictk/bin_table.hpp"
 #include "hictk/cooler/cooler.hpp"
 #include "hictk/cooler/pixel_selector.hpp"
@@ -30,13 +31,15 @@ class PixelSelector {
   using PixelSelectorVar =
       std::variant<cooler::PixelSelector, hic::PixelSelector, hic::PixelSelectorAll>;
   PixelSelectorVar _sel{cooler::PixelSelector{}};
+  std::shared_ptr<const balancing::Weights> _weights{};
 
  public:
   template <typename N>
   class iterator;
 
   template <typename PixelSelectorT>
-  explicit PixelSelector(PixelSelectorT selector);
+  explicit PixelSelector(PixelSelectorT selector,
+                         std::shared_ptr<const balancing::Weights> weights);
 
   template <typename N>
   [[nodiscard]] auto begin(bool sorted = true) const -> iterator<N>;
@@ -59,6 +62,8 @@ class PixelSelector {
 
   [[nodiscard]] PixelSelector fetch(PixelCoordinates coord1_, PixelCoordinates coord2_) const;
 
+  [[nodiscard]] const balancing::Weights &weights() const noexcept;
+
   template <typename PixelSelectorT>
   [[nodiscard]] constexpr const PixelSelectorT &get() const noexcept;
   template <typename PixelSelectorT>
@@ -159,6 +164,8 @@ class File {
   [[nodiscard]] bool has_normalization(std::string_view normalization) const;
   [[nodiscard]] std::vector<balancing::Method> avail_normalizations() const;
   [[nodiscard]] const balancing::Weights &normalization(std::string_view normalization_) const;
+  [[nodiscard]] std::shared_ptr<const balancing::Weights> normalization_ptr(
+      std::string_view normalization_) const;
 
   template <typename FileT>
   [[nodiscard]] constexpr const FileT &get() const noexcept;

src/libhictk/file/include/hictk/impl/file_impl.hpp

@@ -33,7 +33,11 @@
 namespace hictk {
 
 template <typename PixelSelectorT>
-inline PixelSelector::PixelSelector(PixelSelectorT selector) : _sel(std::move(selector)) {}
+inline PixelSelector::PixelSelector(PixelSelectorT selector,
+                                    std::shared_ptr<const balancing::Weights> weights)
+    : _sel(std::move(selector)), _weights(std::move(weights)) {
+  assert(_weights);
+}
 
 template <typename N>
 inline auto PixelSelector::begin([[maybe_unused]] bool sorted) const -> iterator<N> {
@@ -114,14 +118,21 @@ inline PixelSelector PixelSelector::fetch(PixelCoordinates coord1_,
       [&](const auto& sel) -> PixelSelector {
         using T = remove_cvref_t<decltype(sel)>;
         if constexpr (std::is_same_v<hic::PixelSelectorAll, T>) {
-          throw std::runtime_error("TODO");  // TODO fixme
+          throw std::runtime_error(
+              "calling fetch() on a PixelSelector instance set up to fetch genome-wide matrices "
+              "from .hic files is not supported");
         } else {
-          return PixelSelector{sel.fetch(coord1_, coord2_)};
+          return PixelSelector{sel.fetch(coord1_, coord2_), _weights};
         }
       },
       _sel);
 }
 
+inline const balancing::Weights& PixelSelector::weights() const noexcept {
+  assert(_weights);
+  return *_weights;
+}
+
 template <typename PixelSelectorT>
 constexpr const PixelSelectorT& PixelSelector::get() const noexcept {
   return std::get<PixelSelectorT>(_sel);
@@ -217,6 +228,9 @@ inline File::File(std::string uri, std::uint32_t resolution, hic::MatrixType typ
                   hic::MatrixUnit unit) {
   const auto [path, grp] = cooler::parse_cooler_uri(uri);
   if (hic::utils::is_hic_file(path)) {
+    if (resolution == 0) {
+      throw std::runtime_error("resolution cannot be 0 when opening .hic files.");
+    }
     *this = File(hic::File(path, resolution, type, unit));
     return;
   }
@@ -294,7 +308,11 @@ inline std::uint64_t File::nchroms() const {
 }
 
 inline PixelSelector File::fetch(const balancing::Method& normalization) const {
-  return std::visit([&](const auto& fp) { return PixelSelector{fp.fetch(normalization)}; }, _fp);
+  return std::visit(
+      [&](const auto& fp) {
+        return PixelSelector{fp.fetch(normalization), fp.normalization_ptr(normalization)};
+      },
+      _fp);
 }
 
 inline PixelSelector File::fetch(std::string_view range, const balancing::Method& normalization,
@@ -312,7 +330,8 @@ inline PixelSelector File::fetch(std::string_view range1, std::string_view range
                                  hictk::File::QUERY_TYPE query_type) const {
   return std::visit(
       [&](const auto& fp) {
-        return PixelSelector{fp.fetch(range1, range2, normalization, query_type)};
+        return PixelSelector{fp.fetch(range1, range2, normalization, query_type),
+                             fp.normalization_ptr(normalization)};
       },
       _fp);
 }
@@ -324,7 +343,8 @@ inline PixelSelector File::fetch(std::string_view chrom1_name, std::uint32_t sta
   return std::visit(
       [&](const auto& fp) {
         return PixelSelector{
-            fp.fetch(chrom1_name, start1, end1, chrom2_name, start2, end2, normalization)};
+            fp.fetch(chrom1_name, start1, end1, chrom2_name, start2, end2, normalization),
+            fp.normalization_ptr(normalization)};
       },
       _fp);
 }

src/libhictk/hic/include/hictk/hic.hpp

@@ -94,6 +94,15 @@ class File {
   [[nodiscard]] const balancing::Weights &normalization(balancing::Method norm) const;
   [[nodiscard]] const balancing::Weights &normalization(std::string_view norm) const;
 
+  [[nodiscard]] std::shared_ptr<const balancing::Weights> normalization_ptr(
+      balancing::Method norm, const Chromosome &chrom) const;
+  [[nodiscard]] std::shared_ptr<const balancing::Weights> normalization_ptr(
+      std::string_view norm, const Chromosome &chrom) const;
+  [[nodiscard]] std::shared_ptr<const balancing::Weights> normalization_ptr(
+      balancing::Method norm) const;
+  [[nodiscard]] std::shared_ptr<const balancing::Weights> normalization_ptr(
+      std::string_view norm) const;
+
   [[nodiscard]] std::vector<double> expected_values(
       const Chromosome &chrom,
       const balancing::Method &normalization_ = balancing::Method::NONE()) const;

src/libhictk/hic/include/hictk/hic/cache.hpp

@@ -92,6 +92,9 @@ class WeightCache {
   [[nodiscard]] auto find_or_emplace(std::uint32_t chrom_id, balancing::Method norm) -> Value;
   [[nodiscard]] auto find_or_emplace(const Chromosome& chrom, balancing::Method norm) -> Value;
 
+  [[nodiscard]] auto at(std::uint32_t chrom_id, balancing::Method norm) const -> Value;
+  [[nodiscard]] auto at(const Chromosome& chrom, balancing::Method norm) const -> Value;
+
   void clear() noexcept;
   [[nodiscard]] std::size_t size() const noexcept;
 };

src/libhictk/hic/include/hictk/hic/impl/hic_file_impl.hpp

@@ -246,23 +246,42 @@ inline PixelSelector File::fetch(std::uint64_t first_bin1, std::uint64_t last_bi
 
 inline const balancing::Weights& File::normalization(balancing::Method norm,
                                                      const Chromosome& chrom) const {
+  const auto w = normalization_ptr(norm, chrom);
+  assert(w);
+  return *w;
+}
+inline const balancing::Weights& File::normalization(std::string_view norm,
+                                                     const Chromosome& chrom) const {
+  const auto w = normalization_ptr(norm, chrom);
+  assert(w);
+  return *w;
+}
+inline const balancing::Weights& File::normalization(balancing::Method norm) const {
+  const auto w = normalization_ptr(norm);
+  assert(w);
+  return *w;
+}
+inline const balancing::Weights& File::normalization(std::string_view norm) const {
+  const auto w = normalization_ptr(norm);
+  assert(w);
+  return *w;
+}
+
+inline std::shared_ptr<const balancing::Weights> File::normalization_ptr(
+    balancing::Method norm, const Chromosome& chrom) const {
   assert(_weight_cache);
   const auto expected_length = (chrom.size() + bins().resolution() - 1) / bins().resolution();
 
   try {
-    HICTK_DISABLE_WARNING_PUSH
-    HICTK_DISABLE_WARNING_DANGLING_REFERENCE
+    // This takes care of populating the weight cache when appropriate
     const auto& weights = fetch(chrom.name(), norm).weights1();
     if (weights.size() != expected_length) {
       throw std::runtime_error(
           fmt::format(FMT_STRING("{} normalization vector for {} appears to be corrupted: "
                                  "expected {} values, found {}"),
                       norm, chrom.name(), expected_length, weights.size()));
     }
-    // Returning this reference is fine, as the reference points to memory held by a shared_ptr
-    // (which should stay alive at least as long as the current File object exists)
-    return weights;
-    HICTK_DISABLE_WARNING_POP
+    return _weight_cache->at(chrom, norm);
   } catch (const std::exception& e) {
     const std::string_view msg{e.what()};
 
@@ -284,24 +303,25 @@ inline const balancing::Weights& File::normalization(balancing::Method norm,
   *weights = balancing::Weights{std::numeric_limits<double>::quiet_NaN(), expected_length,
                                 balancing::Weights::Type::DIVISIVE};
 
-  return *weights;
+  return weights;
 }
 
-inline const balancing::Weights& File::normalization(std::string_view norm,
-                                                     const Chromosome& chrom) const {
-  return normalization(balancing::Method{norm}, chrom);
+inline std::shared_ptr<const balancing::Weights> File::normalization_ptr(
+    std::string_view norm, const Chromosome& chrom) const {
+  return normalization_ptr(balancing::Method{norm}, chrom);
 }
 
-inline const balancing::Weights& File::normalization(balancing::Method norm) const {
+inline std::shared_ptr<const balancing::Weights> File::normalization_ptr(
+    balancing::Method norm) const {
   assert(_weight_cache);
   auto weights = _weight_cache->find_or_emplace(0, norm);
   if (!weights->empty()) {
-    return *weights;
+    return weights;
   }
 
   if (norm == balancing::Method::NONE()) {
     *weights = balancing::Weights{1.0, bins().size(), balancing::Weights::Type::DIVISIVE};
-    return *weights;
+    return weights;
   }
 
   std::vector<double> buff(bins().size(), std::numeric_limits<double>::quiet_NaN());
@@ -317,11 +337,12 @@ inline const balancing::Weights& File::normalization(balancing::Method norm) con
   }
 
   *weights = balancing::Weights{std::move(buff), balancing::Weights::Type::DIVISIVE};
-  return *weights;
+  return weights;
 }
 
-inline const balancing::Weights& File::normalization(std::string_view norm) const {
-  return normalization(balancing::Method{norm});
+inline std::shared_ptr<const balancing::Weights> File::normalization_ptr(
+    std::string_view norm) const {
+  return normalization_ptr(balancing::Method{norm});
 }
 
 inline std::vector<double> File::expected_values(const Chromosome& chrom,

src/libhictk/hic/include/hictk/hic/impl/weight_cache_impl.hpp

@@ -29,6 +29,14 @@ inline auto WeightCache::find_or_emplace(const Chromosome &chrom, balancing::Met
   return find_or_emplace(chrom.id(), norm);
 }
 
+inline auto WeightCache::at(std::uint32_t chrom_id, balancing::Method norm) const -> Value {
+  return _weights.at(std::make_pair(chrom_id, norm));
+}
+
+inline auto WeightCache::at(const Chromosome &chrom, balancing::Method norm) const -> Value {
+  return at(chrom.id(), norm);
+}
+
 inline void WeightCache::clear() noexcept { _weights.clear(); }
 inline std::size_t WeightCache::size() const noexcept { return _weights.size(); }