Partitions.h

Go to the documentation of this file.

#ifndef LARCOREALG_GEOMETRY_PARTITIONS_H
#define LARCOREALG_GEOMETRY_PARTITIONS_H

// LArSoft libraries
#include "larcorealg/Geometry/SimpleGeo.h" // lar::util::simple_geo::Rectangle
#include "larcorealg/CoreUtils/DebugUtils.h" // lar::debug::demangle(), ...

// C/C++ standard libraries
#include <algorithm> // std::upper_bound()
#include <vector>
#include <string>
#include <sstream>
#include <iterator> // std::next(), ...
#include <memory> // std::unique_ptr<>
#include <utility> // std::move(), std::forward()
#include <type_traits> // std::declval(), std::true_type...
#include <cassert>
#include <cstdlib> // std::size_t


namespace geo {
  
  namespace part {
    
    
    //-------------------------------------------------------------------------
    class AreaOwner {
      
        public:
      using Area_t = lar::util::simple_geo::Rectangle<double>;
      
      using AreaRangeMember_t = Area_t::Range_t (Area_t::*);
      
      AreaOwner(Area_t const& area): myArea(area) {}
      
      bool contains(double w, double d) const
        { return area().contains(w, d); }
      
      Area_t const& area() const { return myArea; }
      
      template <typename Stream>
      void dumpArea(Stream&& out) const { std::forward<Stream>(out) << area(); }
      
        private:
      Area_t myArea; 
      
    }; // class AreaOwner
    
    
    namespace details {
      
      template <AreaOwner::AreaRangeMember_t Range>
      struct PartitionSorterByAreaRangeLower;
      
    } // namespace details
    
    //*************************************************************************
    //***  Fundamental classes
    //-------------------------------------------------------------------------

    template <typename Data>
    struct PartitionDataDescriber {
      template <typename Stream>
      PartitionDataDescriber(
        Stream&& out, Data const* data,
        std::string indent = "", std::string firstIndent = ""
        );
    }; // struct PartitionDataDescriber<>
    
    template <typename Stream, typename Data>
    void describePartitionData(
      Stream&& out, Data const* data,
      std::string indent = "", std::string firstIndent = ""
      );
    
    
    //--------------------------------------------------------------------------
    class PartitionBase: public AreaOwner {
      
        public:
      // imported types
      using Area_t = AreaOwner::Area_t;
      using AreaRangeMember_t = AreaOwner::AreaRangeMember_t;
      
      PartitionBase(Area_t const& area): AreaOwner(area) {}
      
        protected:
      
      std::string describeArea
        (std::string indent, std::string firstIndent) const;
      
    }; // class PartitionBase
    
    
    //--------------------------------------------------------------------------
    template <typename Data>
    class Partition: public PartitionBase {
      
        public:
      using Data_t = Data; 
      using Partition_t = Partition<Data>; 
      using Area_t = PartitionBase::Area_t; 
      
      using Subpartitions_t = std::vector<std::unique_ptr<Partition_t const>>;
      
      Partition(Area_t const& area): PartitionBase(area) {}
      
      virtual ~Partition() = default;
      
      virtual Data_t* data() const { return nullptr; }
      
      virtual Data_t* atPoint(double w, double d) const = 0;
      
      std::string describe
        (std::string indent, std::string firstIndent) const
        { return doDescribe(indent, firstIndent); }
      
      std::string describe(std::string indent = "") const
        { return describe(indent, indent); }
      
      template <typename Pred>
      void walk(Pred&& pred) const { walk(this, pred); }
      
      std::size_t nParts() const { return parts().size(); }
      
        protected:
      
      static Subpartitions_t const NoSubparts; 
      
      virtual Subpartitions_t const& parts() const { return NoSubparts; }
      
      virtual std::string doDescribe
        (std::string indent, std::string firstIndent) const
        { return PartitionBase::describeArea(indent, firstIndent); }
      
      template <typename Pred>
      static void walk(Partition_t const* start, Pred&& pred);
      
    }; // class Partition<>

    //**************************************************************************
    //***  Partition class hierarchy
    //--------------------------------------------------------------------------
    template <typename Data>
    class PartitionWithData: public Partition<Data> {
      
        public:
      using Base_t = Partition<Data>; 
      using Partition_t = Partition<Data>; 
      using Data_t = typename Partition_t::Data_t; 
      using Area_t = typename Partition_t::Area_t; 
      
      PartitionWithData(Area_t const& area, Data_t* myData)
        : Base_t(area), myData(myData)
        {}
      
      virtual Data_t* data() const override { return myData; }
      
      virtual Data_t* atPoint(double w, double d) const override
        { return Base_t::contains(w, d)? myData: nullptr; }
      
        private:
      Data_t* myData; 
      
      virtual std::string doDescribe
        (std::string indent, std::string firstIndent) const override;
      
    }; // class PartitionWithData


    //--------------------------------------------------------------------------
    template <typename Data>
    class PartitionElement: public PartitionWithData<Data> {
      
        public:
      using Base_t = PartitionWithData<Data>; 
      
      // Import constructors
      using Base_t::Base_t;
      
    }; // class PartitionElement


    //--------------------------------------------------------------------------
    template <typename Data>
    class PartitionContainer: public PartitionWithData<Data> {
      
        public:
      using Base_t = PartitionWithData<Data>; 
      using Partition_t = Partition<Data>; 
      
      // inherited types
      using Data_t = typename Partition_t::Data_t;
      using Area_t = typename Partition_t::Area_t;
      using Subpartitions_t = typename Partition_t::Subpartitions_t;
      
      virtual Data_t* atPoint(double w, double d) const override;
      
        protected:
      Subpartitions_t myParts; 
      
      std::size_t size() const { return parts().size(); }
      
      virtual Subpartitions_t const& parts() const override { return myParts; }
      
      PartitionContainer(
        Area_t const& area,
        Subpartitions_t&& subpartitions,
        Data_t* defData = nullptr
        )
        : Base_t(area, defData), myParts(std::move(subpartitions))
        {}
      
      virtual Partition_t const* findPart(double w, double d) const = 0;
      
      virtual std::string doDescribe
        (std::string indent, std::string firstIndent) const override;
      
      virtual std::string describeIntro() const;
      
    }; // class PartitionContainer<>
    
    
    //--------------------------------------------------------------------------
    template <typename Data, typename Sorter>
    class SortedPartition: public PartitionContainer<Data> {
      
        public:
      using Base_t = PartitionContainer<Data>; 
      using Partition_t = Partition<Data>; 
      using Sorter_t = Sorter; 
      
      // inherited data types
      using Data_t = typename Partition_t::Data_t;
      using Area_t = typename Partition_t::Area_t;
      using Subpartitions_t = typename Base_t::Subpartitions_t;
      
      SortedPartition(
        Area_t const& area,
        Subpartitions_t&& subpartitions,
        Data_t* defData = nullptr,
        Sorter_t sorter = {}
        )
        : Base_t(area, std::move(subpartitions), defData), sorter(sorter)
        { initParts(); }
      
        protected:
      
      Sorter_t sorter; 
      
      Partition_t const* findPartWithKey(double key) const;
      
      void initParts();
      
    }; // class SortedPartition<>
    
    
    //--------------------------------------------------------------------------
    template<typename Data, PartitionBase::AreaRangeMember_t Range>
    class PartitionSortedByRange: public SortedPartition
      <Data, details::PartitionSorterByAreaRangeLower<Range>>
    {
        public:
      using Base_t = SortedPartition
        <Data, details::PartitionSorterByAreaRangeLower<Range>>;
      using Base_t::Base_t; // import inherited constructors
      
    }; // class PartitionSortedByRange


    //--------------------------------------------------------------------------
    template <typename Data>
    class DepthPartition
      : public PartitionSortedByRange<Data, &PartitionBase::Area_t::depth>
    {
        public:
      using Base_t
        = PartitionSortedByRange<Data, &PartitionBase::Area_t::depth>;
      
      using Base_t::Base_t; // import inherited constructors
      
      virtual typename Base_t::Partition_t const* findPart
        (double /* w */, double d) const override
        { return Base_t::findPartWithKey(d); }
      
        private:
      
      virtual std::string describeIntro() const override;
      
    }; // class DepthPartition
    
    
    //--------------------------------------------------------------------------
    template <typename Data>
    class WidthPartition
      : public PartitionSortedByRange<Data, &PartitionBase::Area_t::width>
    {
        public:
      using Base_t
        = PartitionSortedByRange<Data, &PartitionBase::Area_t::width>;
      
      using Base_t::Base_t; // inherited constructors
      
      virtual typename Base_t::Partition_t const* findPart
        (double w, double /* d */) const override
        { return Base_t::findPartWithKey(w); }
      
        private:
      
      virtual std::string describeIntro() const override;
      
    }; // class WidthPartition
    
    
    //--------------------------------------------------------------------------
    template <typename Data>
    class GridPartition: public PartitionContainer<Data> {
        public:
      using Partition_t = Partition<Data>; 
      using Base_t = PartitionContainer<Data>; 
      
      // Inherited types
      using Data_t = typename Partition_t::Data_t; 
      using Area_t = typename Partition_t::Area_t; 
      using Subpartitions_t = typename Base_t::Subpartitions_t;
      
      GridPartition(
        Area_t const& area,
        Subpartitions_t&& subpartitions,
        unsigned int nDepthPartitions, unsigned int nWidthPartitions, 
        Data_t* defData = nullptr
        );
      
      GridPartition(
        Area_t const& area,
        Subpartitions_t&& subpartitions,
        unsigned int nDepthPartitions, Data_t* defData = nullptr
        );
      
      
        private:
      std::vector<double> widthSeps; 
      std::vector<double> depthSeps; 
      
      std::size_t nWidthParts() const { return widthSeps.size(); }
      std::size_t nDepthParts() const { return depthSeps.size(); }
      
      auto part(std::size_t iDepth, std::size_t iWidth) ->decltype(auto)
        { return Base_t::parts()[iDepth * nWidthParts() + iWidth]; }
      auto part(std::size_t iDepth, std::size_t iWidth) const ->decltype(auto)
        { return Base_t::parts()[iDepth * nWidthParts() + iWidth]; }
      
      virtual Partition_t const* findPart(double w, double d) const override;
      
      std::vector<double> computeWidthSeps
        (unsigned int nD, unsigned int nW) const;
      std::vector<double> computeDepthSeps
        (unsigned int nD, unsigned int nW) const;
      
      virtual std::string doDescribe
        (std::string indent, std::string firstIndent) const override;
      
      template<
        PartitionBase::AreaRangeMember_t Range,
        typename BeginIter, typename EndIter
        >
      static std::vector<double> detectSeparators(
        BeginIter b, EndIter e,
        std::size_t const nGroups,
        std::size_t const startDelta,
        std::size_t const stride
        );
      
    }; // class GridPartition<>
    
    
    //--------------------------------------------------------------------------
    
  } // namespace part
} // namespace geo


//******************************************************************************
//***  inline implementation
//******************************************************************************
inline std::string geo::part::PartitionBase::describeArea
  (std::string, std::string firstIndent) const
{
  std::ostringstream sstr;
  sstr << firstIndent << "partition covers ";
  dumpArea(sstr);
  return sstr.str();
} // geo::part::PartitionBase::describeArea()


//******************************************************************************
//***  template implementation
//******************************************************************************

namespace geo {
  namespace part {
    namespace details {
      
      //------------------------------------------------------------------------
      //---  some metaprogramming utilities
      //------------------------------------------------------------------------
      template <typename T, typename = void>
      struct is_partition_type: public std::false_type {};
      
      template <typename Part>
      struct is_partition_type
        <
          Part,
          std::enable_if_t
            <std::is_base_of<PartitionBase, std::decay_t<Part>>::value>
        >
        : public std::true_type
      {};
      
      template <typename T>
      constexpr bool is_partition_type_v = is_partition_type<T>();
      
      
      //------------------------------------------------------------------------
      template <typename, typename = void>
      struct is_partition_ptr: public std::false_type {};
      
      template <typename PartPtr>
      struct is_partition_ptr
        <
          PartPtr,
          std::enable_if_t
            <is_partition_type_v<decltype(*std::declval<PartPtr>())>>
        >
        : public std::true_type
      {};
      
      template <typename T>
      constexpr bool is_partition_ptr_v = is_partition_ptr<T>();
      
      
      //------------------------------------------------------------------------
      template <typename, typename = void>
      struct is_partition_ptr_iterator: public std::false_type {};
      
      template <typename Iter>
      struct is_partition_ptr_iterator
        <
          Iter,
          std::enable_if_t
            <is_partition_ptr_v<std::decay_t<typename Iter::value_type>>>
        >
        : public std::true_type
      {};
      
      
      //------------------------------------------------------------------------
      //---  Sorting objects
      //------------------------------------------------------------------------
      template <AreaOwner::AreaRangeMember_t Range>
      struct RangeLowerBoundExtractor {
        static constexpr auto range = Range; 
        using Area_t = AreaOwner::Area_t;
        
        double operator() (double lower) const { return lower; }
        double operator() (Area_t::Range_t const& r) const 
          { return (*this)(r.lower); }
        double operator() (Area_t const& area) const 
          { return (*this)(area.*range); }
        double operator() (AreaOwner const& area) const 
          { return (*this)(area.area()); }
        double operator() (AreaOwner const* ptr) const 
          { return (*this)(*ptr); }
        
      }; // struct RangeLowerBoundExtractor<>
      
      
      //------------------------------------------------------------------------
      template <PartitionBase::AreaRangeMember_t Range>
      struct PartitionRangeLowerBoundExtractor
        : public RangeLowerBoundExtractor<Range>
      {
        using Base_t = RangeLowerBoundExtractor<Range>;
        using Area_t = typename Base_t::Area_t;
        
        using Partition_t = PartitionBase; 
        
        using Base_t::operator(); // import inherited versions
        
        auto operator()(Partition_t const& part)
          { return Base_t::operator()(part.area()); }
        template <
          typename PartPtr,
          typename = std::enable_if_t<details::is_partition_ptr<PartPtr>::value>
          >
        auto operator()(PartPtr const& part)
          { return operator()(*part); }
        
      }; // struct PartitionRangeLowerBoundExtractor<>
      
      
      //------------------------------------------------------------------------
      template <PartitionBase::AreaRangeMember_t Range>
      struct PartitionSorterByAreaRangeLower {
        using Sorter_t = PartitionSorterByAreaRangeLower<Range>;
        
        using KeyExtractor_t = PartitionRangeLowerBoundExtractor<Range>;
        
        static constexpr auto range = Range;
        
        template <typename T>
        static auto key(T const& obj) { return KeyExtractor_t()(obj); }
        
        using Key_t = decltype(Sorter_t::key(std::declval<PartitionBase>()));
        
        static Key_t key(Key_t k) { return k; } // shortcut
        static bool sortKey(Key_t a, Key_t b) { return a < b; }
        
        template <typename A, typename B>
        bool operator() (A const& a, B const& b) const
          { return sortKey(key(a), key(b)); }
        
      }; // struct PartitionSorterByAreaRangeLower
      
      
      //------------------------------------------------------------------------
      
    } // namespace details
  } // namespace part
} // namespace geo


//------------------------------------------------------------------------------
//---  partition data description
//---

template <typename Data>
template <typename Stream>
geo::part::PartitionDataDescriber<Data>::PartitionDataDescriber(
  Stream&& out, Data const* data,
  std::string indent /* = "" */, std::string firstIndent /* = "" */
  )
{
  out << firstIndent;
  std::string typeName = lar::debug::demangle<Data>();
  if (data) {
    out << typeName << "[" << ((void*) data) << "]";
  }
  else { 
    out << "no '" << typeName << "' data";
  }
} // geo::part::PartitionDataDescriber::PartitionDataDescriber()


template <typename Stream, typename Data>
void geo::part::describePartitionData(
  Stream&& out, Data const* data,
  std::string indent /* = "" */, std::string firstIndent /* = "" */
  )
{
  geo::part::PartitionDataDescriber<Data>
    (std::forward<Stream>(out), data, indent, firstIndent);
} // geo::part::describePartitionData()


//------------------------------------------------------------------------------
//--- geo::part::Partition
//---
template <typename Data>
typename geo::part::Partition<Data>::Subpartitions_t const
geo::part::Partition<Data>::NoSubparts;

//------------------------------------------------------------------------------
template <typename Data>
template <typename Pred>
void geo::part::Partition<Data>::walk(Partition_t const* start, Pred&& pred) {
  if (!start) return;
  pred(*start);
  
  // recursive implementation
  for (auto const& subPart: start->parts())
    subPart->walk(std::forward<Pred>(pred));
  
} // geo::part::Partition<Data>::walk()


//------------------------------------------------------------------------------
//---  geo::part::PartitionWithData
//---
template <typename Data>
std::string geo::part::PartitionWithData<Data>::doDescribe
  (std::string indent, std::string firstIndent) const
{
  std::string msg = Base_t::doDescribe(indent, firstIndent);
  if (data()) {
    std::ostringstream sstr; 
    sstr << ": ";
    describePartitionData(sstr, data(), indent);
    msg += sstr.str();
  }
  else {
    msg += " (no data)";
  }
  return msg;
} // geo::part::PartitionWithData<Data>::doDescribe()


//------------------------------------------------------------------------------
//--- geo::part::PartitionContainer
//--- 
//------------------------------------------------------------------------------
template <typename Data>
auto geo::part::PartitionContainer<Data>::atPoint(double w, double d) const
   -> Data_t*
{
  if (!Base_t::contains(w, d)) return nullptr; // not our point at all
  // it's ours; see if it belongs to a subpart
  auto part = findPart(w, d);
  return part? part->atPoint(w, d): Base_t::data();
} // geo::part::PartitionContainer<Data>::atPoint()


//------------------------------------------------------------------------------
template <typename Data>
std::string geo::part::PartitionContainer<Data>::doDescribe
  (std::string indent, std::string firstIndent) const
{
  std::string msg = firstIndent + describeIntro();
  if (Base_t::data()) {
    std::ostringstream sstr; 
    sstr << ", and ";
    describePartitionData(sstr, Base_t::data(), indent, "");
    msg += sstr.str();
  }
  
  for (auto const& part: parts()) {
    msg += "\n" + indent + " * ";
    msg += part->describe(indent + "  ", "");
  }
  
  return msg;
} // geo::part::PartitionContainer<Data>::doDescribe()


//------------------------------------------------------------------------------
template <typename Data>
std::string geo::part::PartitionContainer<Data>::describeIntro() const {
  return std::to_string(parts().size()) + " subpartitions";
} // geo::part::PartitionContainer<Data>::describeIntro()


//------------------------------------------------------------------------------
//---  geo::part::SortedPartition
//---
template <typename Data, typename Sorter>
auto geo::part::SortedPartition<Data, Sorter>::findPartWithKey(double key) const
  -> Partition_t const*
{
  auto pbegin = Base_t::parts().cbegin();
  auto iPart = std::upper_bound(pbegin, Base_t::parts().cend(), key, sorter);
  return (iPart == pbegin)? nullptr: (--iPart)->get();
} // geo::part::SortedPartition<>::findPartWithKey()


//------------------------------------------------------------------------------
template <typename Data, typename Sorter>
void geo::part::SortedPartition<Data, Sorter>::initParts() {
  /*
   * Initialization tasks:
   * - ensure that the parts are sorted by increasing depth
   * 
   */
  std::sort(Base_t::myParts.begin(), Base_t::myParts.end(), sorter);
  
} // geo::part::SortedPartition<>::initParts()


//------------------------------------------------------------------------------
//---  geo::part::DepthPartition
//---
template <typename Data>
std::string geo::part::DepthPartition<Data>::describeIntro() const {
  std::ostringstream sstr;
  sstr
    << Base_t::size() << " partitions along depth covering " << Base_t::area();
  return sstr.str();
} // geo::part::DepthPartition<Data>::describeIntro()


//------------------------------------------------------------------------------
//---  geo::part::WidthPartition
//---
template <typename Data>
std::string geo::part::WidthPartition<Data>::describeIntro() const {
  std::ostringstream sstr;
  sstr
    << Base_t::size() << " partitions along width covering " << Base_t::area();
  return sstr.str();
} // geo::part::WidthPartition<Data>::describeIntro()


//------------------------------------------------------------------------------
//---  geo::part::GridPartition
//---
template <typename Data>
geo::part::GridPartition<Data>::GridPartition(
  Area_t const& area,
  Subpartitions_t&& subpartitions,
    unsigned int nDepthPartitions, unsigned int nWidthPartitions, 
  Data_t* defData /* = nullptr */
  )
  : Base_t(area, std::move(subpartitions), defData)
  , widthSeps(computeWidthSeps(nDepthPartitions, nWidthPartitions))
  , depthSeps(computeDepthSeps(nDepthPartitions, nWidthPartitions))
{
  assert(nWidthPartitions * nDepthPartitions == Base_t::size());
} // geo::part::GridPartition<Data>::GridPartition()


template <typename Data>
geo::part::GridPartition<Data>::GridPartition(
  Area_t const& area,
  Subpartitions_t&& subpartitions,
  unsigned int nDepthPartitions, Data_t* defData /* = nullptr */
  )
  : GridPartition(
    area, std::move(subpartitions), nDepthPartitions,
    (nDepthPartitions? subpartitions.size()/nDepthPartitions: 0),
    defData
    )
  {}


//------------------------------------------------------------------------------
template <typename Data>
std::string geo::part::GridPartition<Data>::doDescribe
  (std::string indent, std::string firstIndent) const
{
  std::ostringstream sstr; 
  sstr << firstIndent << Base_t::describeIntro()
    << " in a (WxD) = " << nWidthParts() << " x " << nDepthParts() << " grid";
  if (Base_t::data()) {
    sstr << ", and ";
    describePartitionData(sstr, Base_t::data(), indent, "");
  }
  for (std::size_t iDepth = 0; iDepth < nDepthParts(); ++iDepth) {
    for (std::size_t iWidth = 0; iWidth < nWidthParts(); ++iWidth) {
      sstr << "\n" << indent << " [" << iDepth << "][" << iWidth << "] "
        << part(iDepth, iWidth)->describe(indent + "  ", "");
    } // for width
  } // for depth
  
  return sstr.str();
} // geo::part::GridPartition<Data>::doDescribe()


//------------------------------------------------------------------------------
template <typename Data>
auto geo::part::GridPartition<Data>::findPart(double w, double d) const
  -> Partition_t const* 
{
  auto const iWidth = std::upper_bound(widthSeps.cbegin(), widthSeps.cend(), w);
  if (iWidth == widthSeps.cbegin()) return nullptr;
  auto const iDepth = std::upper_bound(depthSeps.cbegin(), depthSeps.cend(), d);
  if (iDepth == depthSeps.cbegin()) return nullptr;
  return part(
    std::distance(depthSeps.cbegin(), iDepth) - 1U,
    std::distance(widthSeps.cbegin(), iWidth) - 1U
    ).get();
} // geo::part::GridPartition<Data>::findPart()


//------------------------------------------------------------------------------
template <typename Data>
std::vector<double> geo::part::GridPartition<Data>::computeWidthSeps
  (unsigned int, unsigned int nW) const
{
  return detectSeparators<&Area_t::width>
    (Base_t::parts().cbegin(), Base_t::parts().cend(), nW, 1U, nW);
} // geo::part::GridPartition<Data>::computeWidthSeps()


template <typename Data>
std::vector<double> geo::part::GridPartition<Data>::computeDepthSeps
  (unsigned int nD, unsigned int nW) const
{
  return detectSeparators<&Area_t::depth>
    (Base_t::parts().cbegin(), Base_t::parts().cend(), nD, nW, 1U);
} // geo::part::GridPartition<Data>::computeDepthSeps()


template <typename Data>
template<
  geo::part::AreaOwner::AreaRangeMember_t Range,
  typename BeginIter, typename EndIter
  >
std::vector<double> geo::part::GridPartition<Data>::detectSeparators(
  BeginIter b, EndIter e,
  std::size_t const nGroups,
  std::size_t const startDelta,
  std::size_t const stride
) {
  /*
   * The iterators are better be random access.
   * The range [b,e[ is considered to be a 2D table stored row after row.
   * This function can operate on rows or columns, given the proper arguments.
   * 
   * The full range is split in nGroups "groups" (e.g. rows or columns).
   * Each group g starts at the element g x startDelta.
   * All members of that group are stride elements far from each other.
   * 
   * Be the table of size (d x w).
   * To process data row by row:
   * - the group is a row
   * - the number of groups is the number of rows: nGroups = d
   * - each group will have (d x w)/nGroups = w elements
   * - the start offset is the number of group times the size of it:
   *   startDelta = w
   * - the elements are contiguous: stride = 1
   * 
   * To process data column by column:
   * - the group is a column
   * - the number of groups is the number of columns: nGroups = w
   * - each group will have (d x w)/nGroups = d elements
   * - the start offset matches the number of group: startDelta = 1
   * - the elements are separated by a full row: stride = w
   * 
   */
  // the separator is on the lower bound of selected range of the partition area
  
  static_assert(details::is_partition_ptr_iterator<BeginIter>(),
    "Begin iterator does not point to a pointer to partition type");
  
  details::PartitionRangeLowerBoundExtractor<Range> lowerBound;
  
  std::size_t const nParts = std::distance(b, e);
  std::size_t const nPartsInGroup = nParts / nGroups;
  
  auto const part
    = [b](std::size_t index){ return std::next(b, index)->get(); };
  
  std::vector<double> seps(nGroups);
  for (size_t g = 0; g < nGroups; ++g) {
  
    double& sep = seps[g];
    
    // indices of an element in the previous and next group, respectively
    std::size_t index = g * startDelta;
    sep = lowerBound(part(index));
    
    std::size_t const iend = index + nPartsInGroup * stride;
    while ((index += stride) < iend) {
      double const l = lowerBound(part(index));
      if (sep > l) sep = l;
    } // while (element)
    
  } // while (groups)
  return seps;
} // geo::part::GridPartition<Data>::detectSeparators()


//******************************************************************************

#endif // LARCOREALG_GEOMETRY_PARTITIONS_H