kmer_hash.hpp

Go to the documentation of this file.

// -----------------------------------------------------------------------------------------------------
// Copyright (c) 2006-2020, Knut Reinert & Freie Universität Berlin
// Copyright (c) 2016-2020, Knut Reinert & MPI für molekulare Genetik
// This file may be used, modified and/or redistributed under the terms of the 3-clause BSD-License
// shipped with this file and also available at: https://github.com/seqan/seqan3/blob/master/LICENSE.md
// -----------------------------------------------------------------------------------------------------
 
#pragma once
 
#include <cmath>
 
#include <seqan3/alphabet/concept.hpp>
#include <seqan3/range/hash.hpp>
#include <seqan3/search/kmer_index/shape.hpp>
 
namespace seqan3::detail
{
 
// ---------------------------------------------------------------------------------------------------------------------
// kmer_hash_view class
// ---------------------------------------------------------------------------------------------------------------------
 
template <std::ranges::view urng_t>
class kmer_hash_view : public std::ranges::view_interface<kmer_hash_view<urng_t>>
{
private:
    static_assert(std::ranges::forward_range<urng_t const>, "The kmer_hash_view only works on forward_ranges");
    static_assert(semialphabet<reference_t<urng_t>>, "The reference type of the underlying range must model "
                  "seqan3::semialphabet.");
 
    urng_t urange;
 
    shape shape_;
 
    template <typename rng_t>
    class shape_iterator
    {
    private:
        using it_t = std::ranges::iterator_t<rng_t>;
        using sentinel_t = std::ranges::sentinel_t<rng_t>;
 
    public:
        using difference_type = typename it_t::difference_type;
        using value_type = size_t;
        using pointer = void;
        using reference = value_type;
        using iterator_category = std::input_iterator_tag;
        using iterator_concept = std::conditional_t<std::contiguous_iterator<it_t>,
                                                     typename std::random_access_iterator_tag,
                                                     iterator_tag_t<it_t>>;
 
        constexpr shape_iterator()                                   = default; 
        constexpr shape_iterator(shape_iterator const &)             = default; 
        constexpr shape_iterator(shape_iterator &&)                  = default; 
        constexpr shape_iterator & operator=(shape_iterator const &) = default; 
        constexpr shape_iterator & operator=(shape_iterator &&)      = default; 
        ~shape_iterator()                                            = default; 
 
        shape_iterator(it_t it_start, shape s_) :
            shape_{s_}, text_left{it_start}, text_right{it_start}
        {
            assert(std::ranges::size(shape_) > 0);
 
            roll_factor = std::pow(sigma, std::ranges::size(shape_) - 1);
 
            hash_full();
        }
 
 
        friend bool operator==(shape_iterator const & lhs, sentinel_t const & rhs) noexcept
        {
            return lhs.text_right == rhs;
        }
 
        friend bool operator==(sentinel_t const & lhs, shape_iterator const & rhs) noexcept
        {
            return lhs == rhs.text_right;
        }
 
        friend bool operator==(shape_iterator const & lhs, shape_iterator const & rhs) noexcept
        {
            return std::tie(lhs.text_right, lhs.shape_) == std::tie(rhs.text_right, rhs.shape_);
        }
 
        friend bool operator!=(shape_iterator const & lhs, sentinel_t const & rhs) noexcept
        {
            return !(lhs == rhs);
        }
 
        friend bool operator!=(sentinel_t const & lhs, shape_iterator const & rhs) noexcept
        {
            return !(lhs == rhs);
        }
 
        friend bool operator!=(shape_iterator const & lhs, shape_iterator const & rhs) noexcept
        {
            return !(lhs == rhs);
        }
 
        friend bool operator<(shape_iterator const & lhs, shape_iterator const & rhs) noexcept
        {
            return (lhs.shape_ <= rhs.shape_) && (lhs.text_right < rhs.text_right);
        }
 
        friend bool operator>(shape_iterator const & lhs, shape_iterator const & rhs) noexcept
        {
            return (lhs.shape_ >= rhs.shape_) && (lhs.text_right > rhs.text_right);
        }
 
        friend bool operator<=(shape_iterator const & lhs, shape_iterator const & rhs) noexcept
        {
            return (lhs.shape_ <= rhs.shape_) && (lhs.text_right <= rhs.text_right);
        }
 
        friend bool operator>=(shape_iterator const & lhs, shape_iterator const & rhs) noexcept
        {
            return (lhs.shape_ >= rhs.shape_) && (lhs.text_right >= rhs.text_right);
        }
 
 
        shape_iterator & operator++() noexcept
        {
            hash_forward();
            return *this;
        }
 
        shape_iterator operator++(int) noexcept
        {
            shape_iterator tmp{*this};
            hash_forward();
            return tmp;
        }
 
        shape_iterator & operator--() noexcept
            requires std::bidirectional_iterator<it_t>
        {
            hash_backward();
            return *this;
        }
 
        shape_iterator operator--(int) noexcept
            requires std::bidirectional_iterator<it_t>
        {
            shape_iterator tmp{*this};
            hash_backward();
            return tmp;
        }
 
        shape_iterator & operator+=(difference_type const skip) noexcept
            requires std::random_access_iterator<it_t>
        {
            hash_forward(skip);
            return *this;
        }
 
        shape_iterator operator+(difference_type const skip) const noexcept
            requires std::random_access_iterator<it_t>
        {
            shape_iterator tmp{*this};
            return tmp += skip;
        }
 
        friend shape_iterator operator+(difference_type const skip, shape_iterator const & it) noexcept
            requires std::random_access_iterator<it_t>
        {
            return it + skip;
        }
 
        shape_iterator & operator-=(difference_type const skip) noexcept
            requires std::random_access_iterator<it_t>
        {
            hash_backward(skip);
            return *this;
        }
 
        shape_iterator operator-(difference_type const skip) const noexcept
            requires std::random_access_iterator<it_t>
        {
            shape_iterator tmp{*this};
            return tmp -= skip;
        }
 
        friend shape_iterator operator-(difference_type const skip, shape_iterator const & it) noexcept
            requires std::random_access_iterator<it_t>
        {
            return it - skip;
        }
 
        difference_type operator-(shape_iterator const & lhs) const noexcept
            requires std::random_access_iterator<it_t>
        {
            return static_cast<difference_type>(text_right - lhs.text_right);
        }
 
        friend difference_type operator-(sentinel_t const & lhs, shape_iterator const & rhs) noexcept
            requires std::sized_sentinel_for<sentinel_t, it_t>
        {
            return static_cast<difference_type>(lhs - rhs.text_right);
        }
 
        friend difference_type operator-(shape_iterator const & lhs, sentinel_t const & rhs) noexcept
            requires std::sized_sentinel_for<it_t, sentinel_t>
        {
            return static_cast<difference_type>(lhs.text_right - rhs);
        }
 
        value_type operator[](difference_type const n)
            requires std::random_access_iterator<it_t>
        {
            text_left += n;
            hash_full();
            return operator*();
        }
 
        value_type operator*() const noexcept
        {
            return hash_value + to_rank(*text_right);
        }
 
    private:
        using alphabet_t = value_type_t<it_t>;
 
        static constexpr auto const sigma{alphabet_size<alphabet_t>};
 
        size_t hash_value{0};
 
        size_t roll_factor{0};
 
        shape shape_;
 
        it_t text_left;
 
        it_t text_right;
 
        void hash_forward()
        {
            if (shape_.all())
            {
                hash_roll_forward();
            }
            else
            {
                std::ranges::advance(text_left,  1);
                hash_full();
            }
        }
 
        void hash_forward(difference_type const skip)
            requires std::random_access_iterator<it_t>
        {
            std::ranges::advance(text_left, skip);
            hash_full();
        }
 
        void hash_backward()
            requires std::bidirectional_iterator<it_t>
        {
            if (shape_.all())
            {
                hash_roll_backward();
            }
            else
            {
                std::ranges::advance(text_left,  -1);
                hash_full();
            }
        }
 
        void hash_backward(difference_type const skip)
        {
            std::ranges::advance(text_left, -skip);
            hash_full();
        }
 
        void hash_full()
        {
            text_right = text_left;
            hash_value = 0;
 
            for (size_t i{0}; i < shape_.size() - 1u; ++i)
            {
                hash_value += shape_[i] * to_rank(*text_right);
                hash_value *= shape_[i] ? sigma : 1;
                std::ranges::advance(text_right, 1);
            }
        }
 
        void hash_roll_forward()
        {
            hash_value -= to_rank(*(text_left)) * roll_factor;
            hash_value += to_rank(*(text_right));
            hash_value *= sigma;
 
            std::ranges::advance(text_left,  1);
            std::ranges::advance(text_right, 1);
        }
 
        void hash_roll_backward()
            requires std::bidirectional_iterator<it_t>
        {
            std::ranges::advance(text_left,  -1);
            std::ranges::advance(text_right, -1);
 
            hash_value /= sigma;
            hash_value -= to_rank(*(text_right));
            hash_value += to_rank(*(text_left)) * roll_factor;
        }
    };
 
public:
    kmer_hash_view()                                       = default; 
    kmer_hash_view(kmer_hash_view const & rhs)             = default; 
    kmer_hash_view(kmer_hash_view && rhs)                  = default; 
    kmer_hash_view & operator=(kmer_hash_view const & rhs) = default; 
    kmer_hash_view & operator=(kmer_hash_view && rhs)      = default; 
    ~kmer_hash_view()                                      = default; 
 
    kmer_hash_view(urng_t urange_, shape const & s_) : urange{std::move(urange_)}, shape_{s_}
    {
        if (shape_.size() > (64 / std::log2(alphabet_size<reference_t<urng_t>>)))
        {
            throw std::invalid_argument{"The chosen shape/alphabet combination is not valid. "
                                        "The alphabet or shape size must be reduced."};
        }
    }
 
    template <typename rng_t>
     requires !std::same_as<remove_cvref_t<rng_t>, kmer_hash_view> &&
              std::ranges::viewable_range<rng_t> &&
              std::constructible_from<urng_t, ranges::ref_view<std::remove_reference_t<rng_t>>>
    kmer_hash_view(rng_t && urange_, shape const & s_) :
        urange{std::views::all(std::forward<rng_t>(urange_))}, shape_{s_}
    {
        if (shape_.size() > (64 / std::log2(alphabet_size<reference_t<urng_t>>)))
        {
            throw std::invalid_argument{"The chosen shape/alphabet combination is not valid. "
                                        "The alphabet or shape size must be reduced."};
        }
    }
 
    auto begin() noexcept
    {
        return shape_iterator<urng_t>{std::ranges::begin(urange), shape_};
    }
 
    auto begin() const noexcept
        requires const_iterable_range<urng_t>
    {
        return shape_iterator<urng_t const>{std::ranges::begin(urange), shape_};
    }
 
    auto cbegin() const noexcept
        requires const_iterable_range<urng_t>
    {
        return begin();
    }
 
    auto end() noexcept
    {
        return std::ranges::end(urange);
    }
 
    auto end() const noexcept
        requires const_iterable_range<urng_t>
    {
        return std::ranges::end(urange);
    }
 
    auto cend() const noexcept
        requires const_iterable_range<urng_t>
    {
        return end();
    }
};
 
template <std::ranges::viewable_range rng_t>
kmer_hash_view(rng_t &&, shape const & shape_) -> kmer_hash_view<std::ranges::all_view<rng_t>>;
 
// ---------------------------------------------------------------------------------------------------------------------
// kmer_hash_fn (adaptor definition)
// ---------------------------------------------------------------------------------------------------------------------
 
struct kmer_hash_fn
{
    constexpr auto operator()(shape const & shape_) const
    {
        return adaptor_from_functor{*this, shape_};
    }
 
    template <std::ranges::range urng_t>
    constexpr auto operator()(urng_t && urange, shape const & shape_) const
    {
        static_assert(std::ranges::viewable_range<urng_t>,
            "The range parameter to views::kmer_hash cannot be a temporary of a non-view range.");
        static_assert(std::ranges::forward_range<urng_t>,
            "The range parameter to views::kmer_hash must model std::ranges::forward_range.");
        static_assert(semialphabet<reference_t<urng_t>>,
            "The range parameter to views::kmer_hash must be over elements of seqan3::semialphabet.");
 
        return kmer_hash_view{std::forward<urng_t>(urange), shape_};
    }
};
 
} // namespace seqan3::detail
 
namespace seqan3::views
{
 
inline constexpr auto kmer_hash = detail::kmer_hash_fn{};
 
 
} // namespace seqan3::views