Line data    Source code

       1             : // SPDX-License-Identifier: BSD-3-Clause
       2             : // Copyright (c) 2023 Scipp contributors (https://github.com/scipp)
       3             : /// @file
       4             : /// @author Simon Heybrock
       5             : #include "scipp/core/element/bin.h"
       6             : #include "scipp/core/element/map_to_bins.h"
       7             : 
       8             : #include "scipp/variable/cumulative.h"
       9             : #include "scipp/variable/reduction.h"
      10             : #include "scipp/variable/shape.h"
      11             : #include "scipp/variable/subspan_view.h"
      12             : #include "scipp/variable/transform.h"
      13             : #include "scipp/variable/util.h"
      14             : 
      15             : #include "bin_detail.h"
      16             : 
      17             : namespace scipp::dataset::bin_detail {
      18             : 
      19             : /// Implementation detail of dataset::bin
      20       30512 : void map_to_bins(Variable &out, const Variable &var, const Variable &offsets,
      21             :                  const Variable &indices) {
      22       30512 :   transform_in_place(out, offsets, var, indices, core::element::bin, "bin");
      23       30512 : }
      24             : 
      25        6181 : Variable make_range(const scipp::index begin, const scipp::index end,
      26             :                     const scipp::index stride, const Dim dim) {
      27       12362 :   return cumsum(broadcast(stride * units::none, {dim, (end - begin) / stride}),
      28       12362 :                 dim, CumSumMode::Exclusive);
      29             : }
      30             : 
      31        9459 : void update_indices_by_binning(Variable &indices, const Variable &key,
      32             :                                const Variable &edges, const bool linspace) {
      33        9459 :   const auto dim = edges.dims().inner();
      34        9459 :   if (!indices.dims().includes(key.dims()))
      35           5 :     throw except::BinEdgeError(
      36          10 :         "Requested binning in dimension '" + to_string(dim) +
      37             :         "' but input contains a bin-edge coordinate with no corresponding "
      38             :         "event-coordinate. Provide an event coordinate or convert the "
      39          10 :         "bin-edge coordinate to a non-edge coordinate.");
      40             : 
      41        9454 :   Variable con_edges;
      42        9454 :   Variable edge_view;
      43             : 
      44        9454 :   if (is_bins(edges)) {
      45         721 :     edge_view = as_subspan_view(edges);
      46             :   } else {
      47        8733 :     con_edges = scipp::variable::as_contiguous(edges, dim);
      48        8733 :     edge_view = subspan_view(con_edges.as_const(), dim);
      49             :   }
      50             : 
      51        9454 :   if (linspace) {
      52       13832 :     variable::transform_in_place(
      53       13832 :         indices, key, edge_view.as_const(),
      54             :         core::element::update_indices_by_binning_linspace,
      55             :         "scipp.bin.update_indices_by_binning_linspace");
      56             :   } else {
      57        5076 :     variable::transform_in_place(
      58        5076 :         indices, key, edge_view.as_const(),
      59             :         core::element::update_indices_by_binning_sorted_edges,
      60             :         "scipp.bin.update_indices_by_binning_sorted_edges");
      61             :   }
      62        9454 : }
      63             : 
      64             : namespace {
      65             : template <class Index>
      66         265 : Variable groups_to_map(const Variable &var, const Dim dim) {
      67         265 :   return variable::transform(subspan_view(var, dim),
      68             :                              core::element::groups_to_map<Index>,
      69         265 :                              "scipp.bin.groups_to_map");
      70             : }
      71             : } // namespace
      72             : 
      73        4908 : void update_indices_by_grouping(Variable &indices, const Variable &key,
      74             :                                 const Variable &groups) {
      75        4908 :   const auto dim = groups.dims().inner();
      76        4908 :   const auto con_groups = scipp::variable::as_contiguous(groups, dim);
      77             : 
      78        9813 :   if ((con_groups.dtype() == dtype<int32_t> ||
      79        9791 :        con_groups.dtype() == dtype<int64_t>) &&
      80        4886 :       con_groups.dims().volume() != 0
      81             :       // We can avoid expensive lookups in std::unordered_map if the groups are
      82             :       // contiguous, by simple subtraction of an offset. This is especially
      83             :       // important when the number of target groups is large since the map
      84             :       // lookup would result in frequent cache misses.
      85       14724 :       && isarange(con_groups, con_groups.dim()).value<bool>()) {
      86             :     const auto ngroup = makeVariable<scipp::index>(
      87        4643 :         Values{con_groups.dims().volume()}, units::none);
      88        4643 :     const auto offset = con_groups.slice({con_groups.dim(), 0});
      89        4643 :     variable::transform_in_place(
      90             :         indices, key, ngroup, offset,
      91             :         core::element::update_indices_by_grouping_contiguous,
      92             :         "scipp.bin.update_indices_by_grouping_contiguous");
      93        4643 :     return;
      94        4643 :   }
      95             : 
      96         265 :   const auto map = (indices.dtype() == dtype<int64_t>)
      97             :                        ? groups_to_map<int64_t>(con_groups, dim)
      98         265 :                        : groups_to_map<int32_t>(con_groups, dim);
      99         265 :   variable::transform_in_place(indices, key, map,
     100             :                                core::element::update_indices_by_grouping,
     101             :                                "scipp.bin.update_indices_by_grouping");
     102        4908 : }
     103             : 
     104          38 : void update_indices_from_existing(Variable &indices, const Dim dim) {
     105          38 :   const scipp::index nbin = indices.dims()[dim];
     106          38 :   const auto index = make_range(0, nbin, 1, dim);
     107          38 :   variable::transform_in_place(indices, index, nbin * units::none,
     108             :                                core::element::update_indices_from_existing,
     109             :                                "scipp.bin.update_indices_from_existing");
     110          38 : }
     111             : 
     112             : /// `sub_bin` is a binned variable with sub-bin indices: new bins within bins
     113        8477 : Variable bin_sizes(const Variable &sub_bin, const Variable &offset,
     114             :                    const Variable &nbin) {
     115             :   return variable::transform(
     116       16954 :       as_subspan_view(sub_bin), offset, nbin, core::element::count_indices,
     117       25431 :       "scipp.bin.bin_sizes"); // transform bins, not bin element
     118             : }
     119             : 
     120             : } // namespace scipp::dataset::bin_detail