minimiser.hpp

Go to the documentation of this file.

// SPDX-FileCopyrightText: 2006-2025 Knut Reinert & Freie Universität Berlin
// SPDX-FileCopyrightText: 2016-2025 Knut Reinert & MPI für molekulare Genetik
// SPDX-License-Identifier: BSD-3-Clause
 
#pragma once
 
#include <algorithm>
#include <deque>
 
#include <seqan3/core/detail/empty_type.hpp>
#include <seqan3/core/range/detail/adaptor_from_functor.hpp>
#include <seqan3/core/range/type_traits.hpp>
#include <seqan3/utility/range/concept.hpp>
#include <seqan3/utility/type_traits/lazy_conditional.hpp>
 
namespace seqan3::detail
{
// ---------------------------------------------------------------------------------------------------------------------
// minimiser_view class
// ---------------------------------------------------------------------------------------------------------------------
 
template <std::ranges::view urng1_t, std::ranges::view urng2_t = std::ranges::empty_view<seqan3::detail::empty_type>>
class minimiser_view : public std::ranges::view_interface<minimiser_view<urng1_t, urng2_t>>
{
private:
    static_assert(std::ranges::forward_range<urng1_t>, "The minimiser_view only works on forward_ranges.");
    static_assert(std::ranges::forward_range<urng2_t>, "The minimiser_view only works on forward_ranges.");
    static_assert(std::totally_ordered<std::ranges::range_reference_t<urng1_t>>,
                  "The reference type of the underlying range must model std::totally_ordered.");
 
    using default_urng2_t = std::ranges::empty_view<seqan3::detail::empty_type>;
 
    static constexpr bool second_range_is_given = !std::same_as<urng2_t, default_urng2_t>;
 
    static_assert(!second_range_is_given
                      || std::totally_ordered_with<std::ranges::range_reference_t<urng1_t>,
                                                   std::ranges::range_reference_t<urng2_t>>,
                  "The reference types of the underlying ranges must model std::totally_ordered_with.");
 
    static constexpr bool const_iterable =
        seqan3::const_iterable_range<urng1_t> && seqan3::const_iterable_range<urng2_t>;
 
    urng1_t urange1{};
    urng2_t urange2{};
 
    size_t window_size{};
 
    template <bool const_range>
    class basic_iterator;
 
    using sentinel = std::default_sentinel_t;
 
public:
    minimiser_view()
        requires std::default_initializable<urng1_t> && std::default_initializable<urng2_t>
    = default;                                                        
    minimiser_view(minimiser_view const & rhs) = default;             
    minimiser_view(minimiser_view && rhs) = default;                  
    minimiser_view & operator=(minimiser_view const & rhs) = default; 
    minimiser_view & operator=(minimiser_view && rhs) = default;      
    ~minimiser_view() = default;                                      
 
    explicit minimiser_view(urng1_t urange1, size_t const window_size) :
        minimiser_view{std::move(urange1), default_urng2_t{}, window_size}
    {}
 
    template <typename other_urng1_t>
        requires (std::ranges::viewable_range<other_urng1_t>
                  && std::constructible_from<urng1_t, std::ranges::ref_view<std::remove_reference_t<other_urng1_t>>>)
    explicit minimiser_view(other_urng1_t && urange1, size_t const window_size) :
        urange1{std::views::all(std::forward<other_urng1_t>(urange1))},
        urange2{default_urng2_t{}},
        window_size{window_size}
    {}
 
    explicit minimiser_view(urng1_t urange1, urng2_t urange2, size_t const window_size) :
        urange1{std::move(urange1)},
        urange2{std::move(urange2)},
        window_size{window_size}
    {
        if constexpr (second_range_is_given)
        {
            if (std::ranges::distance(urange1) != std::ranges::distance(urange2))
                throw std::invalid_argument{"The two ranges do not have the same size."};
        }
    }
 
    template <typename other_urng1_t, typename other_urng2_t>
        requires (std::ranges::viewable_range<other_urng1_t>
                  && std::constructible_from<urng1_t, std::views::all_t<other_urng1_t>>
                  && std::ranges::viewable_range<other_urng2_t>
                  && std::constructible_from<urng2_t, std::views::all_t<other_urng2_t>>)
    explicit minimiser_view(other_urng1_t && urange1, other_urng2_t && urange2, size_t const window_size) :
        urange1{std::views::all(std::forward<other_urng1_t>(urange1))},
        urange2{std::views::all(std::forward<other_urng2_t>(urange2))},
        window_size{window_size}
    {
        if constexpr (second_range_is_given)
        {
            if (std::ranges::distance(urange1) != std::ranges::distance(urange2))
                throw std::invalid_argument{"The two ranges do not have the same size."};
        }
    }
 
    basic_iterator<false> begin()
    {
        return {std::ranges::begin(urange1), std::ranges::end(urange1), std::ranges::begin(urange2), window_size};
    }
 
    basic_iterator<true> begin() const
        requires const_iterable
    {
        return {std::ranges::begin(urange1), std::ranges::end(urange1), std::ranges::begin(urange2), window_size};
    }
 
    sentinel end() const
    {
        return {};
    }
};
 
template <std::ranges::view urng1_t, std::ranges::view urng2_t>
template <bool const_range>
class minimiser_view<urng1_t, urng2_t>::basic_iterator
{
private:
    using urng1_sentinel_t = maybe_const_sentinel_t<const_range, urng1_t>;
    using urng1_iterator_t = maybe_const_iterator_t<const_range, urng1_t>;
    using urng2_iterator_t = maybe_const_iterator_t<const_range, urng2_t>;
 
    template <bool>
    friend class basic_iterator;
 
public:
    using difference_type = std::ranges::range_difference_t<urng1_t>;
    using value_type = std::ranges::range_value_t<urng1_t>;
    using pointer = void;
    using reference = value_type;
    using iterator_category = std::forward_iterator_tag;
    using iterator_concept = iterator_category;
 
    basic_iterator() = default;                                   
    basic_iterator(basic_iterator const &) = default;             
    basic_iterator(basic_iterator &&) = default;                  
    basic_iterator & operator=(basic_iterator const &) = default; 
    basic_iterator & operator=(basic_iterator &&) = default;      
    ~basic_iterator() = default;                                  
 
    basic_iterator(basic_iterator<!const_range> const & it)
        requires const_range
        :
        minimiser_value{std::move(it.minimiser_value)},
        urng1_iterator{std::move(it.urng1_iterator)},
        urng1_sentinel{std::move(it.urng1_sentinel)},
        urng2_iterator{std::move(it.urng2_iterator)},
        window_values{std::move(it.window_values)}
    {}
 
    basic_iterator(urng1_iterator_t urng1_iterator,
                   urng1_sentinel_t urng1_sentinel,
                   urng2_iterator_t urng2_iterator,
                   size_t window_size) :
        urng1_iterator{std::move(urng1_iterator)},
        urng1_sentinel{std::move(urng1_sentinel)},
        urng2_iterator{std::move(urng2_iterator)}
    {
        size_t size = std::ranges::distance(urng1_iterator, urng1_sentinel);
        window_size = std::min<size_t>(window_size, size);
 
        window_first(window_size);
    }
 
 
    friend bool operator==(basic_iterator const & lhs, basic_iterator const & rhs)
    {
        return (lhs.urng1_iterator == rhs.urng1_iterator) && (rhs.urng2_iterator == rhs.urng2_iterator)
            && (lhs.window_values.size() == rhs.window_values.size());
    }
 
    friend bool operator!=(basic_iterator const & lhs, basic_iterator const & rhs)
    {
        return !(lhs == rhs);
    }
 
    friend bool operator==(basic_iterator const & lhs, sentinel const &)
    {
        return lhs.urng1_iterator == lhs.urng1_sentinel;
    }
 
    friend bool operator==(sentinel const & lhs, basic_iterator const & rhs)
    {
        return rhs == lhs;
    }
 
    friend bool operator!=(sentinel const & lhs, basic_iterator const & rhs)
    {
        return !(lhs == rhs);
    }
 
    friend bool operator!=(basic_iterator const & lhs, sentinel const & rhs)
    {
        return !(lhs == rhs);
    }
 
    basic_iterator & operator++() noexcept
    {
        next_unique_minimiser();
        return *this;
    }
 
    basic_iterator operator++(int) noexcept
    {
        basic_iterator tmp{*this};
        next_unique_minimiser();
        return tmp;
    }
 
    value_type operator*() const noexcept
    {
        return minimiser_value;
    }
 
    constexpr urng1_iterator_t const & base() const & noexcept
    {
        return urng1_iterator;
    }
 
    constexpr urng1_iterator_t base() &&
    {
        return std::move(urng1_iterator);
    }
 
private:
    value_type minimiser_value{};
 
    size_t minimiser_position_offset{};
 
    urng1_iterator_t urng1_iterator{};
    urng1_sentinel_t urng1_sentinel{};
    urng2_iterator_t urng2_iterator{};
 
    std::deque<value_type> window_values{};
 
    void next_unique_minimiser()
    {
        while (!next_minimiser())
        {}
    }
 
    auto window_value() const
    {
        if constexpr (!second_range_is_given)
            return *urng1_iterator;
        else
            return std::min(*urng1_iterator, *urng2_iterator);
    }
 
    void advance_window()
    {
        ++urng1_iterator;
        if constexpr (second_range_is_given)
            ++urng2_iterator;
    }
 
    void window_first(size_t const window_size)
    {
        if (window_size == 0u)
            return;
 
        for (size_t i = 0u; i < window_size - 1u; ++i)
        {
            window_values.push_back(window_value());
            advance_window();
        }
        window_values.push_back(window_value());
        auto minimiser_it = std::ranges::min_element(window_values, std::less_equal<value_type>{});
        minimiser_value = *minimiser_it;
        minimiser_position_offset = std::distance(std::begin(window_values), minimiser_it);
    }
 
    bool next_minimiser()
    {
        advance_window();
        if (urng1_iterator == urng1_sentinel)
            return true;
 
        value_type const new_value = window_value();
 
        window_values.pop_front();
        window_values.push_back(new_value);
 
        if (minimiser_position_offset == 0)
        {
            auto minimiser_it = std::ranges::min_element(window_values, std::less_equal<value_type>{});
            minimiser_value = *minimiser_it;
            minimiser_position_offset = std::distance(std::begin(window_values), minimiser_it);
            return true;
        }
 
        if (new_value < minimiser_value)
        {
            minimiser_value = new_value;
            minimiser_position_offset = window_values.size() - 1;
            return true;
        }
 
        --minimiser_position_offset;
        return false;
    }
};
 
template <std::ranges::viewable_range rng1_t>
minimiser_view(rng1_t &&, size_t const window_size) -> minimiser_view<std::views::all_t<rng1_t>>;
 
template <std::ranges::viewable_range rng1_t, std::ranges::viewable_range rng2_t>
minimiser_view(rng1_t &&,
               rng2_t &&,
               size_t const window_size) -> minimiser_view<std::views::all_t<rng1_t>, std::views::all_t<rng2_t>>;
 
// ---------------------------------------------------------------------------------------------------------------------
// minimiser_fn (adaptor definition)
// ---------------------------------------------------------------------------------------------------------------------
 
struct minimiser_fn
{
    constexpr auto operator()(size_t const window_size) const
    {
        return adaptor_from_functor{*this, window_size};
    }
 
    template <std::ranges::range urng1_t>
    constexpr auto operator()(urng1_t && urange1, size_t const window_size) const
    {
        static_assert(std::ranges::viewable_range<urng1_t>,
                      "The range parameter to views::minimiser cannot be a temporary of a non-view range.");
        static_assert(std::ranges::forward_range<urng1_t>,
                      "The range parameter to views::minimiser must model std::ranges::forward_range.");
 
        if (window_size == 1) // Would just return urange1 without any changes
            throw std::invalid_argument{"The chosen window_size is not valid. "
                                        "Please choose a value greater than 1 or use two ranges."};
 
        return minimiser_view{urange1, window_size};
    }
};
 
} // namespace seqan3::detail
 
namespace seqan3::views
{
inline constexpr auto minimiser = detail::minimiser_fn{};
 
} // namespace seqan3::views