edit_distance_score_matrix_full.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 <bitset>
 
#include <seqan3/alignment/pairwise/edit_distance_fwd.hpp>
#include <seqan3/core/bit_manipulation.hpp>
 
namespace seqan3::detail
{
 
template <typename word_t, typename score_t, bool is_semi_global, bool use_max_errors>
class edit_distance_score_matrix_full
{
public:
    template <std::ranges::viewable_range database_t,
              std::ranges::viewable_range query_t,
              typename align_config_t,
              typename edit_traits>
    friend class edit_distance_unbanded;
 
     edit_distance_score_matrix_full() = default;                                                    
     edit_distance_score_matrix_full(edit_distance_score_matrix_full const &) = default;             
     edit_distance_score_matrix_full(edit_distance_score_matrix_full &&) = default;                  
     edit_distance_score_matrix_full & operator=(edit_distance_score_matrix_full const &) = default; 
     edit_distance_score_matrix_full & operator=(edit_distance_score_matrix_full &&) = default;      
     ~edit_distance_score_matrix_full() = default;                                                   
 
protected:
    template <typename derived_t, typename edit_traits>
    friend class edit_distance_unbanded_score_matrix_policy;
 
    edit_distance_score_matrix_full(size_t const rows_size)
        : rows_size{rows_size}, columns{}
    {}
 
public:
    using word_type = word_t;
 
    static constexpr auto word_size = sizeof_bits<word_type>;
 
    using score_type = score_t;
 
    using value_type = std::conditional_t<use_max_errors, std::optional<score_type>, score_type>;
 
    using reference = value_type;
 
    using size_type = size_t;
 
    static constexpr std::optional<score_type> inf = std::nullopt;
 
    void reserve(size_t const new_capacity)
    {
        columns.reserve(new_capacity);
    }
 
    template <typename score_type>
    static size_t max_rows(word_type const score_mask, unsigned const last_block,
                           score_type const score, score_type const max_errors) noexcept
    {
        size_t const offset = score_mask == 0u ? 0u : most_significant_bit_set(score_mask) + 1u;
        size_t const active_row = word_size * last_block + offset;
        return active_row + (score <= max_errors);
    }
 
    static score_type score_delta_of_word(word_type const & vp, word_type const & vn) noexcept
    {
        score_type const p = std::bitset<word_size>{vp}.count();
        score_type const n = std::bitset<word_size>{vn}.count();
        return p - n;
    }
 
public:
    reference at(matrix_coordinate const & coordinate) const noexcept
    {
        size_t col = coordinate.col;
        size_t row = coordinate.row;
 
        assert(row < rows());
        assert(col < cols());
 
        column_type const & column = columns[col];
        if constexpr(use_max_errors)
            if (!(row < column.max_rows))
                return inf;
 
        score_type score = is_semi_global ? 0u : static_cast<score_type>(col);
 
        size_t current_row = 1u;
        size_t word_idx = 0u;
 
        for (; current_row + word_size <= row; ++word_idx, current_row += word_size)
            score += score_delta_of_word(column.vp[word_idx], column.vn[word_idx]);
 
        if (row >= current_row)
        {
            word_type const mask = (1u << (row - current_row + 1u)) - 1u;
            score += score_delta_of_word(column.vp[word_idx] & mask, column.vn[word_idx] & mask);
        }
 
        return -score;
    }
 
    size_t rows() const noexcept
    {
        return rows_size;
    }
 
    size_t cols() const noexcept
    {
        return columns.size();
    }
 
protected:
    struct max_errors_state
    {
        size_t max_rows;
    };
 
    struct score_matrix_state
    {
        std::vector<word_type> vp{};
        std::vector<word_type> vn{};
    };
 
    struct column_type :
        enable_state_t<true, score_matrix_state>,
        enable_state_t<use_max_errors, max_errors_state>
    {};
 
    void add_column(std::vector<word_type> vp, std::vector<word_type> vn)
        requires (!use_max_errors)
    {
        column_type column{};
        column.vp = std::move(vp);
        column.vn = std::move(vn);
 
        columns.push_back(std::move(column));
    }
 
    void add_column(std::vector<word_type> vp, std::vector<word_type> vn, size_t const max_rows)
        requires use_max_errors
    {
        column_type column{};
        column.vp = std::move(vp);
        column.vn = std::move(vn);
        column.max_rows = max_rows;
 
        columns.push_back(std::move(column));
    }
 
private:
    size_t rows_size{};
    std::vector<column_type> columns{};
};
 
} // namespace seqan3::detail