SeqAn3  3.2.0
The Modern C++ library for sequence analysis.
alignment_trace_matrix_full_banded.hpp
Go to the documentation of this file.
1 // -----------------------------------------------------------------------------------------------------
2 // Copyright (c) 2006-2022, Knut Reinert & Freie Universität Berlin
3 // Copyright (c) 2016-2022, 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 <iterator>
16 #include <ranges>
17 
24 
25 namespace seqan3::detail
26 {
27 
52 template <typename trace_t, bool coordinate_only = false>
53 class alignment_trace_matrix_full_banded :
54  protected alignment_trace_matrix_base<trace_t>,
55  public alignment_matrix_column_major_range_base<alignment_trace_matrix_full_banded<trace_t, coordinate_only>>
56 {
57 private:
58  static_assert(std::same_as<trace_t, trace_directions> || simd_concept<trace_t>,
59  "Value type must either be a trace_directions object or a simd vector.");
60 
62  using matrix_base_t = alignment_trace_matrix_base<trace_t>;
64  using range_base_t =
65  alignment_matrix_column_major_range_base<alignment_trace_matrix_full_banded<trace_t, coordinate_only>>;
67  friend range_base_t;
68 
69 protected:
70  using typename matrix_base_t::coordinate_type;
71  using typename matrix_base_t::element_type;
72  using typename range_base_t::alignment_column_type;
74  using column_data_view_type =
75  std::conditional_t<coordinate_only,
76  decltype(std::views::iota(coordinate_type{}, coordinate_type{})),
77  decltype(views::zip(std::declval<std::span<element_type>>(),
78  std::declval<std::span<element_type>>(),
79  std::views::iota(coordinate_type{}, coordinate_type{})))>;
80 
81 public:
86  using value_type =
87  alignment_trace_matrix_proxy<coordinate_type,
90  using reference = value_type;
92  using iterator = typename range_base_t::iterator;
94  using sentinel = typename range_base_t::sentinel;
95  using typename matrix_base_t::size_type;
97 
102  constexpr alignment_trace_matrix_full_banded() = default;
104  constexpr alignment_trace_matrix_full_banded(alignment_trace_matrix_full_banded const &) = default;
106  constexpr alignment_trace_matrix_full_banded(alignment_trace_matrix_full_banded &&) = default;
108  constexpr alignment_trace_matrix_full_banded & operator=(alignment_trace_matrix_full_banded const &) = default;
110  constexpr alignment_trace_matrix_full_banded & operator=(alignment_trace_matrix_full_banded &&) = default;
112  ~alignment_trace_matrix_full_banded() = default;
113 
129  template <std::ranges::forward_range first_sequence_t, std::ranges::forward_range second_sequence_t>
130  constexpr alignment_trace_matrix_full_banded(first_sequence_t && first,
131  second_sequence_t && second,
132  align_cfg::band_fixed_size const & band,
133  [[maybe_unused]] trace_t const initial_value = trace_t{})
134  {
135  matrix_base_t::num_cols = static_cast<size_type>(std::ranges::distance(first) + 1);
136  matrix_base_t::num_rows = static_cast<size_type>(std::ranges::distance(second) + 1);
137 
138  band_col_index = std::min<int32_t>(std::max<int32_t>(band.upper_diagonal, 0), matrix_base_t::num_cols - 1);
139  band_row_index =
140  std::min<int32_t>(std::abs(std::min<int32_t>(band.lower_diagonal, 0)), matrix_base_t::num_rows - 1);
141  band_size = band_col_index + band_row_index + 1;
142 
143  // Reserve one more cell to deal with last cell in the banded column which needs only the diagonal and up cell.
144  if constexpr (!coordinate_only)
145  {
146  matrix_base_t::data = typename matrix_base_t::pool_type{number_rows{static_cast<size_type>(band_size)},
147  number_cols{matrix_base_t::num_cols}};
148  matrix_base_t::cache_left.resize(band_size + 1, initial_value);
149  }
150  }
152 
154  auto trace_path(matrix_coordinate const & trace_begin)
155  {
156  static_assert(!coordinate_only, "Requested trace but storing the trace was disabled!");
157 
158  using matrix_iter_t = std::ranges::iterator_t<typename matrix_base_t::pool_type>;
159  using trace_iterator_t = trace_iterator_banded<matrix_iter_t>;
160  using path_t = std::ranges::subrange<trace_iterator_t, std::default_sentinel_t>;
161 
162  if (trace_begin.row >= static_cast<size_t>(band_size) || trace_begin.col >= matrix_base_t::num_cols)
163  throw std::invalid_argument{"The given coordinate exceeds the trace matrix size."};
164 
165  return path_t{trace_iterator_t{matrix_base_t::data.begin() + matrix_offset{trace_begin},
166  column_index_type{band_col_index}},
167  std::default_sentinel};
168  }
169 
171  int32_t band_col_index{};
173  int32_t band_row_index{};
175  int32_t band_size{};
176 
177 private:
179  constexpr alignment_column_type initialise_column(size_type const column_index) noexcept
180  {
181  int32_t slice_begin = std::max<int32_t>(0, band_col_index - column_index);
182  int32_t row_end_index = column_index - band_col_index + band_size;
183  int32_t slice_end = band_size - std::max<int32_t>(row_end_index - matrix_base_t::num_rows, 0);
184 
185  assert(row_end_index >= 0);
186  assert(slice_begin >= 0);
187  assert(slice_end > 0);
188  assert(slice_begin < slice_end);
189 
190  coordinate_type row_begin{column_index_type{column_index}, row_index_type{static_cast<size_type>(slice_begin)}};
191  coordinate_type row_end{column_index_type{column_index}, row_index_type{static_cast<size_type>(slice_end)}};
192  if constexpr (coordinate_only)
193  {
194  return alignment_column_type{
195  *this,
196  column_data_view_type{std::views::iota(std::move(row_begin), std::move(row_end))}};
197  }
198  else
199  {
200  matrix_coordinate band_begin{row_index_type{static_cast<size_type>(slice_begin)},
201  column_index_type{column_index}};
202  size_type slice_size = slice_end - slice_begin;
203  // We need to jump to the offset.
204  auto col =
205  views::zip(std::span<element_type>{std::addressof(matrix_base_t::data[band_begin]), slice_size},
206  std::span<element_type>{std::addressof(matrix_base_t::cache_left[slice_begin]), slice_size},
207  std::views::iota(std::move(row_begin), std::move(row_end)));
208  return alignment_column_type{*this, column_data_view_type{std::move(col)}};
209  }
210  }
211 
213  template <std::random_access_iterator iter_t>
214  constexpr value_type make_proxy(iter_t host_iter) noexcept
215  {
216  if constexpr (coordinate_only)
217  {
218  return {*host_iter, std::ignore, std::ignore, std::ignore, std::ignore};
219  }
220  else
221  {
222  return {
223  std::get<2>(*host_iter), // the current coordinate.
224  std::get<0>(*host_iter), // the current cell.
225  std::get<1>(*(host_iter + 1)), // the last left cell to read from.
226  std::get<1>(*host_iter), // the next left cell to write to.
227  matrix_base_t::cache_up // the last up cell to read/write from/to.
228  };
229  }
230  }
231 };
232 
233 } // namespace seqan3::detail
T addressof(T... args)
Provides seqan3::detail::align_config_band.
Provides seqan3::detail::alignment_matrix_column_major_range_base.
Provides seqan3::detail::alignment_trace_matrix_base.
Provides seqan3::detail::alignment_trace_matrix_proxy.
T declval(T... args)
constexpr auto zip
A view adaptor that takes several views and returns tuple-like values from every i-th element of each...
Definition: zip.hpp:573
matrix_index< size_t > matrix_coordinate
A coordinate type to access an element inside of a two-dimensional matrix.
Definition: matrix_coordinate.hpp:178
The <ranges> header from C++20's standard library.
Provides seqan3::detail::trace_iterator_banded.
Provides seqan3::views::zip.