SeqAn3
The Modern C++ library for sequence analysis.
affine_gap_policy.hpp
Go to the documentation of this file.
1 // -----------------------------------------------------------------------------------------------------
2 // Copyright (c) 2006-2019, Knut Reinert & Freie Universität Berlin
3 // Copyright (c) 2016-2019, Knut Reinert & MPI für molekulare Genetik
4 // This file may be used, modified and/or redistributed under the terms of the 3-clause BSD-License
5 // shipped with this file and also available at: https://github.com/seqan/seqan3/blob/master/LICENSE.md
6 // -----------------------------------------------------------------------------------------------------
7 
13 #pragma once
14 
15 #include <limits>
16 #include <tuple>
17 
20 
21 namespace seqan3::detail
22 {
23 
24 // ----------------------------------------------------------------------------
25 // affine_gap_policy
26 // ----------------------------------------------------------------------------
27 
33 template <typename derived_t, typename cell_t, typename align_local_t = std::false_type>
34 class affine_gap_policy
35 {
36 private:
37 
39  friend derived_t;
40 
44  using score_t = std::tuple_element_t<0, cell_t>;
47 
52  constexpr affine_gap_policy() noexcept = default;
53  constexpr affine_gap_policy(affine_gap_policy const &) noexcept = default;
54  constexpr affine_gap_policy(affine_gap_policy &&) noexcept = default;
55  constexpr affine_gap_policy & operator=(affine_gap_policy const &) noexcept = default;
56  constexpr affine_gap_policy & operator=(affine_gap_policy &&) noexcept = default;
57  ~affine_gap_policy() noexcept = default;
58 
67  template <typename matrix_entry_type, typename cache_t>
68  constexpr void compute_cell(matrix_entry_type && current_cell,
69  cache_t & cache,
70  score_t const score) const noexcept
71  {
72  using std::get;
73 
74  // Unpack the cached variables.
75  auto & [score_entry, coordinate, trace_value] = current_cell;
76  auto & [main_score, hz_score, hz_trace] = score_entry;
77  auto & [prev_cell, gap_open, gap_extend, opt] = cache;
78  auto & [prev_score, vt_score, vt_trace] = prev_cell;
79 
80  //TODO For the local alignment we might be able to use GCC overflow builtin arithmetics, which
81  // allows us to check if overflow/underflow would happen. Not sure, if this helps with the performance though.
82  // (see https://gcc.gnu.org/onlinedocs/gcc/Integer-Overflow-Builtins.html)
83 
84  // Precompute the diagonal score.
85  score_t tmp = prev_score + score;
86 
87  // Compute where the max comes from.
88  if constexpr (decays_to_ignore_v<decltype(trace_value)>) // Don't compute a traceback
89  {
90  tmp = (tmp < vt_score) ? vt_score : tmp;
91  tmp = (tmp < hz_score) ? hz_score : tmp;
92  if constexpr (align_local_t::value)
93  tmp = (tmp < 0) ? 0 : tmp;
94  }
95  else // Compute any traceback
96  {
97  tmp = (tmp < vt_score) ? (trace_value = vt_trace, vt_score)
98  : (trace_value = trace_directions::diagonal | vt_trace, tmp);
99  tmp = (tmp < hz_score) ? (trace_value = hz_trace, hz_score)
100  : (trace_value |= hz_trace, tmp);
101  if constexpr (align_local_t::value)
102  tmp = (tmp < 0) ? (trace_value = trace_directions::none, 0) : tmp;
103  }
104  // Cache the current main score for the next diagonal computation and update the current score.
105  prev_score = main_score;
106  main_score = tmp;
107  // Check if this was the optimum. Possibly a noop.
108  static_cast<derived_t const &>(*this).check_score(
109  alignment_optimum{tmp, static_cast<alignment_coordinate>(coordinate)}, opt);
110 
111  // Prepare horizontal and vertical score for next column.
112  tmp += gap_open;
113  vt_score += gap_extend;
114  hz_score += gap_extend;
115 
116  if constexpr (decays_to_ignore_v<decltype(trace_value)>)
117  {
118  vt_score = (vt_score < tmp) ? tmp : vt_score;
119  hz_score = (hz_score < tmp) ? tmp : hz_score;
120  }
121  else
122  {
123  vt_score = (vt_score < tmp) ? (vt_trace = trace_directions::up_open, tmp)
124  : (vt_trace = trace_directions::up, vt_score);
125  hz_score = (hz_score < tmp) ? (hz_trace = trace_directions::left_open, tmp)
126  : (hz_trace = trace_directions::left, hz_score);
127  }
128  }
129 
135  template <typename gap_scheme_t>
136  constexpr auto make_cache(gap_scheme_t && scheme) const noexcept
137  {
138  alignment_optimum opt{std::numeric_limits<score_t>::lowest(),
139  alignment_coordinate{column_index_type{0u}, row_index_type{0u}}};
140  return std::tuple{cell_t{},
141  static_cast<score_t>(scheme.get_gap_open_score() + scheme.get_gap_score()),
142  static_cast<score_t>(scheme.get_gap_score()),
143  std::move(opt)};
144  }
145 };
146 
147 } // namespace seqan3::detail
T lowest(T... args)
Provides the declaration of seqan3::detail::trace_directions.
Definition: aligned_sequence_concept.hpp:35
auto const get
A view calling std::get on each element in a range.
Definition: get.hpp:66
Provides seqan3::detail::alignment_optimum.