diff --git a/CMakeLists.txt b/CMakeLists.txt index 94c0204..86a5e1e 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -104,6 +104,7 @@ target_link_libraries(fastq_bgzf_file PUBLIC add_library(bam_header OBJECT src/bam_header.cpp) +add_library(bamrec OBJECT src/bamrec.cpp) add_library(bam_file OBJECT src/bam_file.cpp) target_link_libraries(bam_file PUBLIC HTSLIB::HTSLIB) @@ -210,6 +211,7 @@ target_link_libraries(falco PUBLIC fastq_gz_file fastq_bgzf_file bam_file + bamrec sam_file samrec bam_header diff --git a/src/bam_file.cpp b/src/bam_file.cpp index 8c7a32c..cd152bc 100644 --- a/src/bam_file.cpp +++ b/src/bam_file.cpp @@ -7,6 +7,7 @@ #include "falco_word.hpp" #include "task_queue.hpp" +#include // IWYU pragma: keep #include #include @@ -40,7 +41,7 @@ estimate_n_reads_bam(const std::string &filename) throw std::runtime_error("failed to identify file format: " + filename); // NOLINTNEXTLINE(cppcoreguidelines-pro-type-union-access) - const auto &fp = format->format == bam ? f->fp.bgzf->fp : f->fp.hfile; + const auto fp = format->format == bam ? f->fp.bgzf->fp : f->fp.hfile; const auto pos_after_header = htell(fp); std::unique_ptr rec(bam_init1(), &bam_destroy1); std::uint64_t n_reads{}; @@ -53,11 +54,12 @@ estimate_n_reads_bam(const std::string &filename) throw std::system_error(std::make_error_code(std::errc(errno)), "error reading bam record from: " + filename); const auto pos_after_reads = htell(fp); - const auto n_compressed_bytes = pos_after_reads - pos_after_header; + const auto n_bytes = pos_after_reads - pos_after_header; const auto filesize = std::filesystem::file_size(filename); - const auto estimate = static_cast( - as_frac(n_reads * (filesize - pos_after_header), n_compressed_bytes)); - return {estimate, total_read_len / n_reads, filesize}; + const auto n_reads_estimate = static_cast( + as_frac(n_reads * (filesize - pos_after_header), n_bytes)); + const auto read_len_estimate = total_read_len / n_reads; + return {n_reads_estimate, read_len_estimate, filesize}; } [[nodiscard]] auto @@ -138,10 +140,10 @@ partition(auto itr, // std::atomic_int32_t &n_tasks) { // ADS: this isn't working as desired: the end position of each part should be // the first record end past the 'end_itr' below unless end_itr == end - const auto dist = std::distance(itr, end); + auto dist = std::distance(itr, end); const auto chunk_size = (dist + n_chunks - 1) / n_chunks; while (itr != end) { - const auto dist = std::distance(itr, end); + dist = std::distance(itr, end); auto end_itr = itr + (dist < chunk_size ? dist : chunk_size); // ADS: find_end_pos doesn't find end pos of a record, but of a range, so // includes multiple records @@ -164,13 +166,14 @@ bam_file::get_chunks(const std::int64_t n_chunks, // // input buffer has been inflated and will provide data for analysis std::swap(input_buffer, output_buffer); std::swap(output_last, input_last); - auto itr = std::data(output_buffer); - const auto end = itr + output_last; // NOLINT(*-pointer-arithmetic) - if (bh) // if we are still parsing the header + const auto beg = std::cbegin(output_buffer); + const auto end = beg + output_last; + auto itr = beg; + if (bh) // if we are still parsing the header itr = bh.update(itr, end); - if (std::distance(itr, end) > 0) + if (itr != end) itr = partition(itr, end, n_chunks, file_id, tq, n_tasks); - output_cursor = std::distance(std::data(output_buffer), itr); + output_cursor = std::distance(beg, itr); } auto diff --git a/src/bam_header.cpp b/src/bam_header.cpp index 1455b7f..5f7bd3e 100644 --- a/src/bam_header.cpp +++ b/src/bam_header.cpp @@ -6,8 +6,8 @@ #include [[nodiscard]] auto -bam_header::update(bam_header::iterator itr, - const bam_header::iterator end) -> bam_header::iterator { +bam_header::update(const_iterator itr, + const const_iterator end) -> const_iterator { static constexpr auto msg = "incorrect BAM magic identified: {} at {}"; const auto update_u32 = [](auto &val, const auto inc, const auto the_byte) { // NOLINTNEXTLINE(*-avoid-magic-numbers) diff --git a/src/bam_header.hpp b/src/bam_header.hpp index 49305dc..8f4068f 100644 --- a/src/bam_header.hpp +++ b/src/bam_header.hpp @@ -7,10 +7,11 @@ #include #include +#include struct bam_header { - using iterator = char *; - using const_iterator = const char *; + using iterator = std::vector::iterator; + using const_iterator = std::vector::const_iterator; static constexpr auto magic = "BAM\1"; static constexpr auto magic_size = 4; @@ -30,7 +31,7 @@ struct bam_header { } [[nodiscard]] auto - update(iterator itr, const iterator end) -> iterator; + update(const_iterator itr, const const_iterator end) -> const_iterator; [[nodiscard]] auto ref_incomplete() const -> bool { diff --git a/src/bamrec.cpp b/src/bamrec.cpp new file mode 100644 index 0000000..287d169 --- /dev/null +++ b/src/bamrec.cpp @@ -0,0 +1,121 @@ +// SPDX-License-Identifier: MIT; Copyright 2026 Andrew D Smith + +#include "bamrec.hpp" + +#include +#include +#include +#include +#include +#include +#include + +[[nodiscard]] auto +bamrec::to_string() const -> std::string { + const auto buffer_s = std::span(std::cbegin(buffer), std::cend(buffer)); + const auto name_itr = std::cbegin(buffer_s); + const auto seq_itr = name_itr + name_len; + const auto qual_itr = seq_itr + seq_len; + std::string qual_fixed(seq_len, '\0'); + std::transform(qual_itr, qual_itr + seq_len, std::begin(qual_fixed), + [](const auto c) { return c + quality_score_offset; }); + return std::format("@{}\n{}\n+\n{}", // + std::string(name_itr, seq_itr), // + std::string(seq_itr, qual_itr), // + qual_fixed); +} + +#ifdef seq_nt16_str +#undef seq_nt16_str +#endif + +#ifdef bam_seqi +#undef bam_seqi +#endif + +template +static inline constexpr OutputIt +assign_sequence_revcomp(BidirIt first, auto last, OutputIt d_first) { + constexpr auto complem = [](const auto x) { + return "TNGNNNCNNNNNNNNNNNNA"[x - 'A']; + }; + constexpr auto seq_nt16_str = "=ACMGRSVTWYHKDBN"; + constexpr auto bam_seqi = [](const auto s, const auto i) -> int { + constexpr auto low_nibble_on = 0xf; + return s[i >> 1] >> ((~i & 1) << 2) & low_nibble_on; + }; + for (auto j = last; j != 0; ++d_first) + // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic) + *d_first = complem(seq_nt16_str[bam_seqi(first, --j)]); + return d_first; +} + +template +static inline constexpr OutputIt +assign_sequence(BidirIt first, auto last, OutputIt d_first) { + constexpr auto seq_nt16_str = "=ACMGRSVTWYHKDBN"; + constexpr auto bam_seqi = [](const auto s, const auto i) -> int { + constexpr auto low_nibble_on = 0xf; + return s[i >> 1] >> ((~i & 1) << 2) & low_nibble_on; + }; + for (auto j = 0U; j != last; ++j) + // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic) + *d_first++ = seq_nt16_str[bam_seqi(first, j)]; + return d_first; +} + +[[nodiscard]] auto +bamrec::get_next(bamrec::pos_t &itr, const bamrec::pos_t end, + bamrec &rec) -> bool { + if (std::distance(itr, end) < bam_core_t::sz) + return false; + bam_core_t core{}; + std::memcpy(std::addressof(core), std::to_address(itr), sizeof(bam_core_t)); + if (std::distance(itr, end) < core.real_block_size()) + return false; + rec.name_len = core.l_read_name - 1; // we don't need the '\0' + rec.seq_len = core.l_seq; + const auto rec_size = rec.name_len + 2 * rec.seq_len; + if (std::size(rec.buffer) < rec_size) + rec.buffer.resize(rec_size); + auto out_itr = std::begin(rec.buffer); + std::copy_n(itr + bam_core_t::read_name_offset, rec.name_len, out_itr); + out_itr += rec.name_len; // increment data cursor to sequence + auto seq_in = itr + core.seq_offset(); + if (core.bam_is_rev()) + assign_sequence_revcomp(seq_in, core.l_seq, out_itr); + else + assign_sequence(seq_in, core.l_seq, out_itr); + out_itr += rec.seq_len; // increment data cursor to qual + const auto has_qual = (itr[core.qual_offset()] != qual_missing_code); + if (has_qual) { + const auto qual_itr = itr + core.qual_offset(); + const auto qual_end = qual_itr + rec.seq_len; + if (core.bam_is_rev()) + std::reverse_copy(qual_itr, qual_end, out_itr); + else + std::copy(qual_itr, qual_end, out_itr); + } + else + *out_itr = static_cast(qual_missing_code); + itr += core.real_block_size(); + return true; +} + +[[nodiscard]] auto +bamrec::find_end_pos(bamrec::pos_t itr, + const bamrec::pos_t end) -> bamrec::pos_t { + static constexpr std::int64_t record_size_size = sizeof(std::uint32_t); + std::uint32_t record_size{}; + while (itr != end) { + if (std::distance(itr, end) < record_size_size) + return itr; + std::memcpy(std::addressof(record_size), std::to_address(itr), + record_size_size); + record_size += record_size_size; + if (std::distance(itr, end) < record_size) + return itr; + itr += record_size; // only increment on consume of full record + } + return itr; +} diff --git a/src/bamrec.hpp b/src/bamrec.hpp index 4d15c08..a2fd618 100644 --- a/src/bamrec.hpp +++ b/src/bamrec.hpp @@ -3,100 +3,72 @@ #ifndef SRC_BAMREC_HPP_ #define SRC_BAMREC_HPP_ -#include // for BGZF -#include // for htell - -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include +#include #include -#include -#include -#include -#include #include -#include -#include // IWYU print: keep -#include #include #ifdef bam_is_rev #undef bam_is_rev #endif -#ifdef bam_seqi -#undef bam_seqi -#endif - #ifdef BAM_FREVERSE #undef BAM_FREVERSE #endif -class bamrec { -public: - using pos_t = char *; - -private: - struct bamrec_core_t { - static constexpr auto sz = 36; - static constexpr auto read_name_offset = 36; - std::uint32_t block_size{}; // 0 - std::int32_t refID{}; // 4 - std::int32_t pos{}; // 8 - std::uint8_t l_read_name{}; // 12 - std::uint8_t mapq{}; // 13 - std::uint16_t bin{}; // 14 - std::uint16_t n_cigar_op{}; // 16 - std::uint16_t flag{}; // 18 - std::uint32_t l_seq{}; // 20 - std::int32_t next_refID{}; // 24 - std::int32_t next_pos{}; // 28 - std::int32_t tlen{}; // 32 - // 36 +struct bam_core_t { + static constexpr auto sz = 36; + static constexpr auto read_name_offset = 36; + std::uint32_t block_size{}; // 0 + std::int32_t refID{}; // 4 + std::int32_t pos{}; // 8 + std::uint8_t l_read_name{}; // 12 + std::uint8_t mapq{}; // 13 + std::uint16_t bin{}; // 14 + std::uint16_t n_cigar_op{}; // 16 + std::uint16_t flag{}; // 18 + std::uint32_t l_seq{}; // 20 + std::int32_t next_refID{}; // 24 + std::int32_t next_pos{}; // 28 + std::int32_t tlen{}; // 32 + // 36 + + [[nodiscard]] constexpr auto + seq_offset() const -> std::uint32_t { + static constexpr auto cigar_op_size = sizeof(std::uint32_t); + return read_name_offset + l_read_name + cigar_op_size * n_cigar_op; + } - [[nodiscard]] constexpr auto - seq_offset() const -> std::uint32_t { - static constexpr auto cigar_op_size = sizeof(std::uint32_t); - return read_name_offset + l_read_name + cigar_op_size * n_cigar_op; - } + [[nodiscard]] constexpr auto + n_seq_bytes() const -> std::uint32_t { + return (l_seq + 1) / 2; + } - [[nodiscard]] constexpr auto - n_seq_bytes() const -> std::uint32_t { - return (l_seq + 1) / 2; - } + [[nodiscard]] constexpr auto + qual_offset() const -> std::uint32_t { + return seq_offset() + n_seq_bytes(); + } - [[nodiscard]] constexpr auto - qual_offset() const -> std::uint32_t { - return seq_offset() + n_seq_bytes(); - } + [[nodiscard]] constexpr auto + real_block_size() const -> std::uint32_t { + return (sizeof block_size) + block_size; + } - [[nodiscard]] constexpr auto - real_block_size() const -> std::uint32_t { - return (sizeof block_size) + block_size; - } + [[nodiscard]] constexpr auto + bam_is_rev() -> bool { + static constexpr auto BAM_FREVERSE = 16; + return (flag & BAM_FREVERSE) != 0; + } +}; // bam_core_t - [[nodiscard]] constexpr auto - bam_is_rev() -> bool { - static constexpr auto BAM_FREVERSE = 16; - return (flag & BAM_FREVERSE) != 0; - } - }; +class bamrec { +public: + using pos_t = std::vector::const_iterator; +private: static constexpr auto quality_score_offset = 33; - static constexpr std::uint8_t qual_missing_code = 0xff; // from sam.c - static constexpr auto l_read_name_offset = 8; - static constexpr auto n_cigar_op_offset = 12; - static constexpr auto flag_offset = 14; - static constexpr auto l_seq_offset = 16; + static constexpr char qual_missing_code = -1; // 0xff from sam.c - bamrec_core_t core; std::vector buffer; std::uint32_t name_len{}; std::uint32_t seq_len{}; @@ -111,9 +83,7 @@ class bamrec { friend constexpr auto get_qual(const bamrec &); friend constexpr auto get_qual_end(const bamrec &); friend constexpr auto get_qual_size(const bamrec &); - // clang-format on - // clang-format off bamrec() = default; ~bamrec() = default; bamrec(const bamrec &) = delete; @@ -123,34 +93,20 @@ class bamrec { // clang-format on [[nodiscard]] auto - to_string() const { - // NOLINTBEGIN(cppcoreguidelines-pro-bounds-pointer-arithmetic) - auto qual_itr = std::data(buffer) + name_len + seq_len; - std::string qual_fixed(seq_len, '\0'); - std::transform(qual_itr, qual_itr + seq_len, std::begin(qual_fixed), - [](const auto c) { return c + quality_score_offset; }); - return std::format( - "@{}\n{}\n+\n{}", - std::string(std::data(buffer), std::data(buffer) + name_len), - std::string(std::data(buffer) + name_len, - std::data(buffer) + name_len + seq_len), - qual_fixed); - // NOLINTEND(cppcoreguidelines-pro-bounds-pointer-arithmetic) - } + to_string() const -> std::string; operator bool() const { return name_len != 0; } [[nodiscard]] static auto - get_next(auto &itr, const auto end, bamrec &rec) -> bool; + get_next(pos_t &itr, const pos_t end, bamrec &rec) -> bool; [[nodiscard]] static auto - find_end_pos(pos_t itr, const pos_t end) -> bamrec::pos_t; + find_end_pos(pos_t itr, const pos_t end) -> pos_t; }; -// NOLINTBEGIN(cppcoreguidelines-pro-bounds-pointer-arithmetic) [[nodiscard]] inline constexpr auto get_name(const bamrec &rec) { - return std::data(rec.buffer); + return std::cbegin(rec.buffer); } [[nodiscard]] inline constexpr auto @@ -187,112 +143,15 @@ get_qual_end(const bamrec &rec) { get_qual_size(const bamrec &rec) { return get_seq_size(rec); } -// NOLINTEND(cppcoreguidelines-pro-bounds-pointer-arithmetic) - -using bam_chunks_t = std::vector>; struct bam_task_t { bamrec::pos_t beg{}; bamrec::pos_t end{}; }; -#ifdef seq_nt16_str -#undef seq_nt16_str -#endif - -#ifdef bam_seqi -#undef bam_seqi -#endif - -template -static inline constexpr OutputIt -assign_sequence_revcomp(BidirIt first, auto last, OutputIt d_first) { - constexpr auto complem = [](const auto x) { - return "TNGNNNCNNNNNNNNNNNNA"[x - 'A']; - }; - constexpr auto seq_nt16_str = "=ACMGRSVTWYHKDBN"; - constexpr auto bam_seqi = [](const auto s, const auto i) -> int { - constexpr auto low_nibble_on = 0xf; - // NOLINTNEXTLINE(cppcoreguidelines-avoid-magic-numbers) - return s[i >> 1] >> ((~i & 1) << 2) & low_nibble_on; - }; - for (auto j = last; j != 0; ++d_first) - *d_first = complem(seq_nt16_str[bam_seqi(first, --j)]); - return d_first; -} - -template -static inline constexpr OutputIt -assign_sequence(BidirIt first, auto last, OutputIt d_first) { - constexpr auto seq_nt16_str = "=ACMGRSVTWYHKDBN"; - constexpr auto bam_seqi = [](const auto s, const auto i) -> int { - constexpr auto low_nibble_on = 0xf; - return s[i >> 1] >> ((~i & 1) << 2) & low_nibble_on; - }; - for (auto j = 0U; j != last; ++j) - *d_first++ = seq_nt16_str[bam_seqi(first, j)]; - return d_first; -} - [[nodiscard]] inline auto -bamrec::get_next(auto &itr, const auto end, bamrec &rec) -> bool { - bamrec::bamrec_core_t &core = rec.core; - if (std::distance(itr, end) < bamrec::bamrec_core_t::sz) - return false; - std::memcpy(&core, itr, sizeof core); - if (std::distance(itr, end) < core.real_block_size()) - return false; - rec.name_len = (core.l_read_name - 1); // we don't need the '\0' - rec.seq_len = core.l_seq; - const auto rec_size = rec.name_len + 2 * rec.seq_len; - if (std::size(rec.buffer) < rec_size) - rec.buffer.resize(rec_size); - auto data_itr = std::data(rec.buffer); - std::memcpy(data_itr, itr + bamrec_core_t::read_name_offset, rec.name_len); - // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic) - data_itr += rec.name_len; // increment data cursor to sequence - auto seq_in = itr + core.seq_offset(); - if (core.bam_is_rev()) - assign_sequence_revcomp(seq_in, core.l_seq, data_itr); - else - assign_sequence(seq_in, core.l_seq, data_itr); - // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic) - data_itr += rec.seq_len; // increment data cursor to qual - const auto has_qual = (itr[core.qual_offset()] != qual_missing_code); - if (has_qual) { - const auto qual_itr = itr + core.qual_offset(); - if (core.bam_is_rev()) - std::reverse_copy(qual_itr, qual_itr + rec.seq_len, data_itr); - else - // std::memcpy(data_itr, itr + core.qual_offset(), rec.seq_len); - std::copy(qual_itr, qual_itr + rec.seq_len, data_itr); - } - else - // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic) - data_itr[0] = static_cast(qual_missing_code); - itr += core.real_block_size(); - return true; -} - -[[nodiscard]] inline auto -get_next(auto &itr, const auto end, bamrec &rec) -> bool { +get_next(bamrec::pos_t &itr, const bamrec::pos_t end, bamrec &rec) -> bool { return bamrec::get_next(itr, end, rec); } -[[nodiscard]] inline auto -bamrec::find_end_pos(pos_t itr, const pos_t end) -> bamrec::pos_t { - static constexpr std::int64_t record_size_size = sizeof(std::uint32_t); - std::uint32_t record_size{}; - while (itr != end) { - if (std::distance(itr, end) < record_size_size) - return itr; - std::memcpy(&record_size, itr, record_size_size); - record_size += record_size_size; - if (std::distance(itr, end) < record_size) - return itr; - itr += record_size; // only increment on consume of full record - } - return itr; -} - #endif // SRC_BAMREC_HPP_ diff --git a/src/bgzf_block.hpp b/src/bgzf_block.hpp index 521b8e2..ef4ff61 100644 --- a/src/bgzf_block.hpp +++ b/src/bgzf_block.hpp @@ -4,7 +4,6 @@ #define SRC_BGZF_BLOCK_HPP_ #include -#include static constexpr auto max_bgzf_block_size = 65536; @@ -21,14 +20,11 @@ struct bgzf_block_t { // clang-format off bgzf_block_t(const bgzf_block_t &src) = delete; auto operator=(const bgzf_block_t &src) -> bgzf_block_t & = delete; - bgzf_block_t() = default; ~bgzf_block_t() = default; bgzf_block_t(bgzf_block_t &&src) noexcept = default; auto operator=(bgzf_block_t &&src) noexcept -> bgzf_block_t & = default; - [[nodiscard]] auto data() const -> const char * { return in_itr; } - [[nodiscard]] auto operator<=>(const bgzf_block_t &other) const = default; operator bool() const { return size > 0; } // clang-format on @@ -36,8 +32,6 @@ struct bgzf_block_t { decompress() -> void; }; -using bgzf_chunks_t = std::vector; - inline auto decompress(bgzf_block_t &x) -> void { // free function x.decompress(); diff --git a/src/bgzf_reader.cpp b/src/bgzf_reader.cpp index 65fe015..1805c5a 100644 --- a/src/bgzf_reader.cpp +++ b/src/bgzf_reader.cpp @@ -1,8 +1,9 @@ // SPDX-License-Identifier: MIT; Copyright 2026 Andrew D Smith #include "bgzf_reader.hpp" -#include "bgzf_block.hpp" // for bgzf_block_t, max_bgzf_block_size +#include "bgzf_block.hpp" +#include #include #include #include @@ -11,31 +12,33 @@ #include #include #include +#include -// NOLINTBEGIN(*-bounds-pointer-arithmetic,*-avoid-magic-numbers,*-type-reinterpret-cast) [[nodiscard]] static inline constexpr auto -get_unaligned_le32(const std::uint8_t *p) -> std::uint32_t { - return (static_cast(p[3]) << 24) | - (static_cast(p[2]) << 16) | - (static_cast(p[1]) << 8) | - (static_cast(p[0]) << 0); +get_unaligned_le32(const auto p) -> std::int32_t { + std::int32_t value{}; + std::memcpy(std::addressof(value), p, sizeof(std::int32_t)); + if constexpr (std::endian::native == std::endian::big) + return std::byteswap(value); + else + return value; } [[nodiscard]] static inline constexpr auto -get_isize(const auto *data, const auto data_size) { +get_isize(const auto data, const auto data_size) { static constexpr decltype(data_size) isize_size = 4; assert(data_size > isize_size); - const auto u_data = reinterpret_cast(data); - const auto u_data_isize = u_data + data_size - isize_size; - return data_size < isize_size - ? 0 - : static_cast(get_unaligned_le32(u_data_isize)); + const auto data_isize = data + data_size - isize_size; + return data_size < isize_size ? 0 : get_unaligned_le32(data_isize); } auto -assign(gzip_header &hdr, auto *data) -> void { - // NOLINTNEXTLINE(*-type-reinterpret-cast) - std::memcpy(reinterpret_cast(&hdr), data, gzip_header_size); +assign(gzip_header &hdr, const auto data) -> void { + // ADS: data from the file takes 18 bytes, the fields of the struct take 18 + // bytes, but the struct occupies 20 due to uint32_t members + // NOLINTNEXTLINE(cppcoreguidelines-pro-type-reinterpret-cast) + std::memcpy(reinterpret_cast(&hdr), std::to_address(data), + gzip_header_size); } bgzf_reader::bgzf_reader(const std::string &filename, @@ -57,15 +60,16 @@ bgzf_reader::read_data() -> bool { const auto unused_in = std::distance(next_in, end_in); std::memcpy(inbuf.get(), next_in, unused_in); next_in = inbuf.get(); + // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic) end_in = next_in + unused_in; const auto avail_in = inbuf_size - unused_in; - const auto r = std::fread(end_in, 1, avail_in, fp.get()); - if (std::ferror(fp.get())) { - const auto errc = std::make_error_code(std::errc(errno)); - throw std::system_error(errc, "failed reading input"); - } - end_in += r; // will usually be end of inbuf - return r > 0; + const auto n_bytes = std::fread(end_in, 1, avail_in, fp.get()); + if (std::ferror(fp.get())) + throw std::system_error(std::make_error_code(std::errc(errno)), + "failed reading input"); + // will usually be end of inbuf + end_in += n_bytes; // NOLINT(cppcoreguidelines-pro-bounds-pointer-arithmetic) + return n_bytes > 0; } [[nodiscard]] inline constexpr auto @@ -81,6 +85,7 @@ bgzf_reader::get_decomp_task(char *out_itr) -> bgzf_block_t { return bgzf_block_t{}; assign(gh, next_in); assert(gh.check_magic()); + // NOLINTBEGIN(cppcoreguidelines-pro-bounds-pointer-arithmetic) next_in += gzip_header_size; const auto body_size = get_gzip_body_size(gh); if (std::distance(next_in, end_in) < body_size && !read_data()) @@ -89,6 +94,6 @@ bgzf_reader::get_decomp_task(char *out_itr) -> bgzf_block_t { bgzf_block_t task(get_isize(next_in, body_size), out_itr, next_out); next_in += body_size; next_out += max_bgzf_block_size; + // NOLINTEND(cppcoreguidelines-pro-bounds-pointer-arithmetic) return task; } -// NOLINTEND(*-bounds-pointer-arithmetic,*-avoid-magic-numbers,*-type-reinterpret-cast) diff --git a/src/bgzf_reader.hpp b/src/bgzf_reader.hpp index da5cbc0..8f3e629 100644 --- a/src/bgzf_reader.hpp +++ b/src/bgzf_reader.hpp @@ -31,6 +31,8 @@ struct gzip_header { std::uint16_t size{}; // 16 [2] // 18 + // ADS: total size is 20 bytes because of alignment and 32-bit values + [[nodiscard]] auto check_magic() const -> bool { return id1 == magic1 && id2 == magic2; // ADS: check b and c also diff --git a/src/falco_utils.cpp b/src/falco_utils.cpp index 2feba81..f150925 100644 --- a/src/falco_utils.cpp +++ b/src/falco_utils.cpp @@ -9,7 +9,6 @@ #include // for std::localtime #include #include // for std::put_time -#include #include #include #include diff --git a/src/fastq_file.cpp b/src/fastq_file.cpp index 4c4face..bf8421b 100644 --- a/src/fastq_file.cpp +++ b/src/fastq_file.cpp @@ -2,15 +2,24 @@ #include "fastq_file.hpp" #include "falco_utils.hpp" +#include "fqrec.hpp" +#include "task_queue.hpp" #include #include #include #include +#include #include +#include +#include +#include #include #include +#include +#include +#include [[nodiscard]] auto estimate_n_reads_fastq(const std::string &filename) @@ -27,34 +36,117 @@ estimate_n_reads_fastq(const std::string &filename) fstat(fd, &buf); const auto filesize = buf.st_size; if (filesize < n_parts) - return {0LU, 0LU, filesize}; + return {{}, {}, filesize}; const auto part_size = filesize < n_parts * max_part_size ? filesize / n_parts : max_part_size; auto total_newlines = 0ul; - auto read_len_est = 0ul; + auto readlen_est = 0ul; for (auto i = 0; i < n_parts; ++i) { const auto offset = (i * part_size) & page_mask; - char *data = static_cast( - mmap(nullptr, part_size, PROT_READ, MAP_PRIVATE, fd, offset)); - if (data == MAP_FAILED) + auto raw = mmap(nullptr, part_size, PROT_READ, MAP_PRIVATE, fd, offset); + if (raw == MAP_FAILED) throw std::system_error(std::make_error_code(std::errc(errno)), "failed to mmap file"); - // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic) - total_newlines += std::ranges::count(data, data + part_size, '\n'); - read_len_est += estimate_read_length_fastq_chunk(data, part_size); - - if (munmap(static_cast(data), part_size)) + const auto data = std::span(static_cast(raw), part_size); + total_newlines += std::ranges::count(data, '\n'); + readlen_est += estimate_read_length_fastq_chunk(std::data(data), part_size); + if (munmap(raw, part_size)) throw std::system_error(std::make_error_code(std::errc(errno)), "failed to unmap memory"); } close(fd); - read_len_est /= n_parts; + readlen_est /= n_parts; const auto n_reads_est = as_frac(total_newlines, fastq_lines_per_read) * as_frac(static_cast(filesize), (part_size * n_parts)); - return {static_cast(n_reads_est), read_len_est, filesize}; + return {static_cast(n_reads_est), readlen_est, filesize}; +} + +static inline auto +mmap_fastq(const int fd, const std::int64_t offset, const std::int64_t length, + auto &data) { + static constexpr auto prot = PROT_READ; + static constexpr auto flags = MAP_PRIVATE; + data = static_cast(mmap(nullptr, length, prot, flags, fd, offset)); + if (data == MAP_FAILED) + throw std::system_error(std::make_error_code(std::errc(errno)), + "failed to mmap file"); +} + +static inline auto +cleanup_mmap_fastq(auto &buffer, std::int64_t &buffer_size) { + if (buffer == nullptr) + return; + munmap(static_cast(buffer), buffer_size); + buffer = nullptr; + buffer_size = 0; +} + +auto +fastq_file::reset() -> void { + cleanup_mmap_fastq(buffer, buffer_size); +} + +auto +fastq_file::load_next() -> void { + // memory mapped data is page aligned but the data we need is not + static const auto page_mask = sysconf(_SC_PAGESIZE) - 1; + std::tie(start_in_file, cursor) = [&] { + const auto pos_in_file = start_in_file + cursor; + return std::tuple(pos_in_file & (~page_mask), pos_in_file & page_mask); + }(); + stop_in_file = std::min(filesize, start_in_file + target_buffer_size); + if (buffer_size > 0) + cleanup_mmap_fastq(buffer, buffer_size); + if (start_in_file < stop_in_file) { // this exist for empty files + buffer_size = stop_in_file - start_in_file; + mmap_fastq(fd, start_in_file, buffer_size, buffer); + } +} + +auto +fastq_file::get_chunks(const std::int64_t n_chunks, const std::int32_t file_id, + task_queue &tq, std::atomic_int32_t &n_tasks) -> void { + static constexpr auto rec_lines = 4; // FASTQ + assert(n_chunks > 0); + const auto buf = std::span(buffer, buffer_size); + // clang-format off + const auto not_read_start = [&](const auto p) { + // ADS: could get confused if '+' lines have full name info + return buf[p] != '@' || (p > 0 && buf[p - 1] != '\n') || + (p > 2 && buf[p - 2] == '+' && buf[p - 3] == '\n'); + }; + const auto fwd_to_read_start = [&](auto pos) { + if (pos == 0) return pos; + while (pos < buffer_size && not_read_start(pos)) ++pos; + return pos; + }; + const auto rev_to_read_start = [&](auto pos) { + while (pos > 0 && (pos == buffer_size || not_read_start(pos))) --pos; + return pos; + }; + // clang-format on + const auto n_bytes_available = buffer_size - cursor; + std::int64_t start_pos = cursor; + const auto [chunk_size, remainder] = std::div(n_bytes_available, n_chunks); + std::vector> chunks(n_chunks); + std::int64_t chunk_end{start_pos}; + for (const auto chunk_idx : std::views::iota(0, n_chunks)) { + const auto chunk_beg = fwd_to_read_start(start_pos); + const auto stop_pos = start_pos + chunk_size + (chunk_idx < remainder); + chunk_end = fwd_to_read_start(stop_pos); + if (chunk_idx + 1 == n_chunks) + if (const auto prev = rev_to_read_start(chunk_end); + std::count(std::cbegin(buf) + prev, std::cend(buf), '\n') < rec_lines) + chunk_end = prev; + ++n_tasks; + // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic) + tq.push(file_id, fq_task_t(buffer + chunk_beg, buffer + chunk_end)); + start_pos = stop_pos; + } + cursor = chunk_end; } diff --git a/src/fastq_file.hpp b/src/fastq_file.hpp index d4ae5a3..624bd06 100644 --- a/src/fastq_file.hpp +++ b/src/fastq_file.hpp @@ -3,70 +3,38 @@ #ifndef SRC_FASTQ_FILE_HPP_ #define SRC_FASTQ_FILE_HPP_ -#include "fqrec.hpp" -#include "task_queue.hpp" - #include -#include #include #include #include -#include #include #include -#include #include #include -#include #include #include #include -#include #include -#include -struct fastq_buffer { - char *data{}; // not necessarily owned - std::int64_t sz{}; // slight redundancy with vars containing classes -}; - -static inline auto -mmap_fastq(const int fd, const std::int64_t start_pos_in_file, - const std::int64_t stop_pos_in_file, fastq_buffer &buf) { - const auto n_bytes = stop_pos_in_file - start_pos_in_file; - char *data = static_cast( - mmap(nullptr, n_bytes, PROT_READ, MAP_PRIVATE, fd, start_pos_in_file)); - if (data == MAP_FAILED) - throw std::system_error(std::make_error_code(std::errc(errno)), - "failed to mmap file"); - buf.data = data; - buf.sz = n_bytes; -} - -static inline auto -cleanup_mmap_fastq(fastq_buffer &buf) { - if (buf.data == nullptr) - return; - [[maybe_unused]] const auto r = munmap(static_cast(buf.data), buf.sz); - buf.data = nullptr; - buf.sz = 0; -} +struct task_queue; struct fastq_file { - static constexpr auto min_buf_size = 64 * 1024; - std::int64_t buf_size{}; + static constexpr auto min_buf_size = 65536; + std::int64_t target_buffer_size{}; std::int64_t filesize{}; - fastq_buffer buf{}; - std::int64_t start_pos_in_file{}; - std::int64_t stop_pos_in_file{}; + char *buffer{}; + std::int64_t buffer_size{}; + std::int64_t start_in_file{}; + std::int64_t stop_in_file{}; std::int64_t cursor{}; int fd{}; - fastq_file(const std::string &filename, const std::int64_t buf_size_arg) : - buf_size{buf_size_arg + min_buf_size}, + fastq_file(const std::string &filename, + const std::int64_t target_buffer_size) : + target_buffer_size{ + std::max(target_buffer_size, static_cast(min_buf_size))}, filesize{static_cast(std::filesystem::file_size(filename))}, - stop_pos_in_file{buf_size_arg + min_buf_size}, // init to use as sentinel fd{open(std::data(filename), O_RDONLY, 0)} { if (fd < 0) throw std::system_error(std::make_error_code(std::errc(errno)), @@ -80,100 +48,34 @@ struct fastq_file { // clang-format on fastq_file(fastq_file &&src) noexcept : - buf_size{src.buf_size}, // - filesize{src.filesize}, // - buf{src.buf}, // - start_pos_in_file{src.start_pos_in_file}, // - stop_pos_in_file{src.stop_pos_in_file}, // - cursor{src.cursor}, // - fd{dup(src.fd)} // <- LOOK + target_buffer_size{src.target_buffer_size}, // + filesize{src.filesize}, // + buffer{src.buffer}, // + buffer_size{src.buffer_size}, // + start_in_file{src.start_in_file}, // + stop_in_file{src.stop_in_file}, // + cursor{src.cursor}, // + fd{dup(src.fd)} // <- LOOK {} auto - reset() -> void { - if (buf.sz > 0) - cleanup_mmap_fastq(buf); - } + reset() -> void; ~fastq_file() { reset(); close(fd); // will always have been opened using a filename } - [[nodiscard]] operator bool() const { - return stop_pos_in_file - start_pos_in_file == buf_size; - } + operator bool() const { return stop_in_file != filesize; } auto - load_next() -> void { - // memory mapped data is page aligned but the data we need is not - static const auto page_mask = sysconf(_SC_PAGESIZE) - 1; - std::tie(start_pos_in_file, cursor) = [&] { - const auto pos_in_file = start_pos_in_file + cursor; - return std::tuple(pos_in_file & (~page_mask), pos_in_file & page_mask); - }(); - stop_pos_in_file = std::min(filesize, start_pos_in_file + buf_size); - if (buf.sz > 0) - cleanup_mmap_fastq(buf); - if (start_pos_in_file < stop_pos_in_file) // this exist for empty files - mmap_fastq(fd, start_pos_in_file, stop_pos_in_file, buf); - } + load_next() -> void; - [[nodiscard]] auto - get_chunks(const std::int64_t n_chunks) -> fq_chunks_t; + auto + get_chunks(const std::int64_t n_chunks, const std::int32_t file_id, + task_queue &tq, std::atomic_int32_t &n_tasks) -> void; }; -[[nodiscard]] static inline auto -get_chunks_fastq_impl(auto &fq, const std::int64_t n_chunks) { - static constexpr auto rec_lines = 4; // FASTQ - // clang-format off - const auto not_read_start = [](const auto s, const auto p) { - // ADS: could get confused if '+' lines have full name info - return s[p] != '@' || (p > 0 && s[p-1] != '\n') || - (p > 2 && s[p-2] == '+' && s[p-3] == '\n'); - }; - const auto fwd_to_read_start = [&](const auto &buf, auto pos) { - if (pos == 0) return pos; - while (pos < buf.sz && not_read_start(buf.data, pos)) ++pos; - return pos; - }; - const auto rev_to_read_start = [&](const auto &buf, auto pos) { - while (pos > 0 && (pos == buf.sz || not_read_start(buf.data, pos))) --pos; - return pos; - }; - // clang-format on - const auto &buf = fq.buf; // non-copyable - const auto n_bytes_available = buf.sz - fq.cursor; - const auto [chunk_size, remainder] = std::div(n_bytes_available, n_chunks); - std::vector> chunks(n_chunks); - std::int64_t start_pos = fq.cursor; - for (const auto chunk_idx : std::views::iota(0, n_chunks)) { - const auto chunk_beg = fwd_to_read_start(buf, start_pos); - const auto stop_pos = start_pos + chunk_size + (chunk_idx < remainder); - const auto chunk_end = fwd_to_read_start(buf, stop_pos); - chunks[chunk_idx] = {chunk_beg, chunk_end}; - start_pos = stop_pos; - } - // make sure final chunk includes only full records - const auto prev_start = rev_to_read_start(buf, chunks.back().second); - // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic) - if (std::count(buf.data + prev_start, buf.data + buf.sz, '\n') < rec_lines) - chunks.back().second = prev_start; - fq.cursor = chunks.back().second; - return chunks; -} - -[[nodiscard]] inline auto -fastq_file::get_chunks(const std::int64_t n_chunks) -> fq_chunks_t { - assert(n_chunks > 0); - const auto add_offset = [d = buf.data](const auto &x) { - // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic) - return std::pair{d + x.first, d + x.second}; - }; - return get_chunks_fastq_impl(*this, n_chunks) | - std::views::transform(add_offset) | std::ranges::to(); -} - [[nodiscard]] auto estimate_n_reads_fastq(const std::string &filename) -> std::tuple; @@ -188,11 +90,7 @@ make_tasks(fastq_file &reads_file, // const auto n_chunks = n_chunks_per_thread * n_threads; n_tasks = 1; // for current task, which makes tasks reads_file.load_next(); - auto chunks = reads_file.get_chunks(n_chunks); - for (const auto &chunk : chunks) { - ++n_tasks; // ADS: increment before submitting - tq.push(file_id, fq_task_t(chunk.first, chunk.second)); - } + reads_file.get_chunks(n_chunks, file_id, tq, n_tasks); } inline auto diff --git a/src/fqrec.hpp b/src/fqrec.hpp index 826f6e5..8bf0357 100644 --- a/src/fqrec.hpp +++ b/src/fqrec.hpp @@ -53,8 +53,6 @@ struct fq_task_t { fqrec::pos_t end{}; }; -using fq_chunks_t = std::vector>; - [[nodiscard]] inline auto get_next(fqrec::pos_t &cursor, const fqrec::pos_t end_itr) -> fqrec { // ADS: need to make sure cursor < end_itr or we will move past diff --git a/src/sam_file.cpp b/src/sam_file.cpp index 76f48a0..f03293d 100644 --- a/src/sam_file.cpp +++ b/src/sam_file.cpp @@ -46,58 +46,52 @@ sam_file::skip_header() -> bool { auto sam_file::shift_output_buffer() -> void { - if (cursor == 0) // shifting at cursor == 0 would do nothing + if (cursor == std::cbegin(buffer)) // shifting here does nothing return; - const auto n_bytes = last - cursor; - std::copy_n(std::cbegin(buffer) + cursor, n_bytes, std::begin(buffer)); - last = n_bytes; - cursor = 0; + last = std::copy(cursor, last, std::begin(buffer)); + cursor = std::begin(buffer); } auto sam_file::get_chunks(const std::int64_t n_chunks, // const std::int32_t file_id, // task_queue &tq, // - std::atomic_int32_t &n_tasks // - ) -> void { + std::atomic_int32_t &n_tasks) -> void { assert(n_chunks > 0); - // NOLINTBEGIN(cppcoreguidelines-pro-bounds-pointer-arithmetic) - - const auto data = std::data(buffer); + const auto beg_itr = std::begin(buffer); + const auto end_itr = last; // clang-format off const auto fwd_to_read_start = [&](auto p) { - if (p == 0) return p; - while (p < last && data[p] != '\n') ++p; - if (p < last) ++p; + if (p == beg_itr) return p; + while (p != end_itr && *p != '\n') ++p; + if (p != end_itr) ++p; return p; }; const auto rev_to_read_start = [&](auto p) { - while (p > 0 && data[p - 1] != '\n') --p; + while (p != beg_itr && *(p - 1) != '\n') --p; return p; }; // clang-format on - const auto n_bytes_available = last; + const auto n_bytes_available = std::distance(beg_itr, end_itr); const auto chunk_size = (n_bytes_available + n_chunks - 1) / n_chunks; assert(n_chunks > 0); - std::int64_t start_pos{}; - std::int64_t chunk_end{}; + auto start_pos = beg_itr; + auto chunk_end = beg_itr; for (const auto chunk_idx : std::views::iota(0, n_chunks)) { const auto chunk_beg = fwd_to_read_start(start_pos); auto stop_pos = start_pos + chunk_size; chunk_end = fwd_to_read_start(stop_pos); - if (chunk_idx == n_chunks - 1) { // final chunk has only full records + if (chunk_idx + 1 == n_chunks) { // final chunk has only full records const auto prev_start = rev_to_read_start(chunk_end); - const auto n_trailing = std::count(data + prev_start, data + last, '\n'); + const auto n_trailing = std::count(prev_start, end_itr, '\n'); if (n_trailing == 0) chunk_end = prev_start; } ++n_tasks; - // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic) - tq.push(file_id, sam_task_t(data + chunk_beg, data + chunk_end)); + tq.push(file_id, sam_task_t(chunk_beg, chunk_end)); start_pos = stop_pos; } cursor = chunk_end; - // NOLINTEND(cppcoreguidelines-pro-bounds-pointer-arithmetic) } auto @@ -105,11 +99,10 @@ sam_file::make_tasks(const std::int64_t n_chunks, // const std::int32_t file_id, // task_queue &tq, // std::atomic_int32_t &n_tasks) -> void { - n_tasks = 1; // for current task, which makes tasks + n_tasks = 1; // for current task, which makes more tasks shift_output_buffer(); - auto data = std::data(buffer); - // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic) - const auto r = std::fread(data + last, 1, std::size(buffer) - last, in.get()); + const auto n_bytes = std::distance(last, std::end(buffer)); + const auto r = std::fread(std::to_address(last), 1, n_bytes, in.get()); if (std::ferror(in.get())) std::system_error(std::make_error_code(std::errc(errno)), "error reading SAM file"); diff --git a/src/sam_file.hpp b/src/sam_file.hpp index f52d5dd..03f547a 100644 --- a/src/sam_file.hpp +++ b/src/sam_file.hpp @@ -4,8 +4,10 @@ #define SRC_SAM_FILE_HPP_ #include +#include #include #include +#include #include #include #include @@ -15,8 +17,8 @@ struct task_queue; class sam_file { static constexpr auto min_buf_size = 64 * 1024; std::vector buffer; - std::int64_t last{}; - std::int64_t cursor{}; + std::vector::iterator cursor; + std::vector::iterator last; std::unique_ptr in; public: @@ -39,6 +41,8 @@ class sam_file { reset() -> void { buffer.clear(); buffer.shrink_to_fit(); + cursor = std::begin(buffer); + last = cursor; } private: diff --git a/src/samrec.cpp b/src/samrec.cpp index 185ee02..7368fde 100644 --- a/src/samrec.cpp +++ b/src/samrec.cpp @@ -6,7 +6,10 @@ #include #include #include +#include #include +#include +#include #include #include @@ -18,6 +21,18 @@ #undef bam_is_rev #endif +[[nodiscard]] auto +samrec::to_string() const -> std::string { + const auto buffer_s = std::span(std::cbegin(buffer), std::cend(buffer)); + const auto name_s = buffer_s.subspan(0, name_len); + const auto seq_s = buffer_s.subspan(name_len, seq_len); + const auto qual_s = buffer_s.subspan(name_len + seq_len, seq_len); + return std::format("@{}\n{}\n+\n{}", // + std::string(std::cbegin(name_s), std::cend(name_s)), + std::string(std::cbegin(seq_s), std::cend(seq_s)), + std::string(std::cbegin(qual_s), std::cend(qual_s))); +} + [[nodiscard]] auto samrec::get_next(samrec::pos_t &cursor, const samrec::pos_t end_itr, samrec &rec) -> bool { @@ -29,13 +44,12 @@ samrec::get_next(samrec::pos_t &cursor, const samrec::pos_t end_itr, static constexpr auto BAM_FREVERSE = 16; return (flag & BAM_FREVERSE) != 0; }; + const auto next_delim = [end_itr](auto &x) { + x = std::find(x, end_itr, '\t'); + }; auto itr = cursor; - const auto next_delim = [end_itr](auto &itr) { - itr = std::find(itr, end_itr, '\t'); - }; - // get the read name const auto name_itr = itr; next_delim(itr); @@ -44,15 +58,13 @@ samrec::get_next(samrec::pos_t &cursor, const samrec::pos_t end_itr, return false; rec.name_len = name_len; - // NOLINTBEGIN(cppcoreguidelines-pro-bounds-pointer-arithmetic) - // get the flag (for revcomp) int flag{}; - auto [flag_end, ec] = std::from_chars(itr, end_itr, flag); + auto [flag_end, ec] = + std::from_chars(std::to_address(itr), std::to_address(end_itr), flag); if (ec != std::errc{}) return false; - // NOLINTNEXTLINE(cppcoreguidelines-pro-type-const-cast) - itr += std::distance(itr, const_cast(flag_end)); + itr += std::distance(std::to_address(itr), flag_end); if (itr++ == end_itr) return false; @@ -80,13 +92,13 @@ samrec::get_next(samrec::pos_t &cursor, const samrec::pos_t end_itr, if (qual_len != rec.seq_len) throw std::runtime_error("quality scores invalid for record: " + - std::string(name_itr, name_len)); + std::string(name_itr, name_itr + name_len)); const auto rec_size = rec.name_len + 2 * rec.seq_len; if (std::ssize(rec.buffer) < rec_size) rec.buffer.resize(rec_size); - auto out_name_itr = std::data(rec.buffer); + auto out_name_itr = std::begin(rec.buffer); auto out_seq_itr = out_name_itr + name_len; auto out_qual_itr = out_seq_itr + seq_len; std::copy_n(name_itr, name_len, out_name_itr); @@ -110,7 +122,5 @@ samrec::get_next(samrec::pos_t &cursor, const samrec::pos_t end_itr, cursor = itr; - // NOLINTEND(cppcoreguidelines-pro-bounds-pointer-arithmetic) - return true; } diff --git a/src/samrec.hpp b/src/samrec.hpp index e69391d..3e96f8d 100644 --- a/src/samrec.hpp +++ b/src/samrec.hpp @@ -4,27 +4,17 @@ #define SRC_SAMREC_HPP_ #include -#include #include #include -#include #include -#ifdef bam_is_rev -#undef bam_is_rev -#endif - -#ifdef BAM_FREVERSE -#undef BAM_FREVERSE -#endif - class samrec { public: - using pos_t = char *; + using pos_t = std::vector::const_iterator; private: - static constexpr auto qual_missing_symbol = '*'; // from SAMv1.pdf - static constexpr std::uint8_t qual_missing_code = 0xff; // from sam.c + static constexpr auto qual_missing_symbol = '*'; // from SAMv1.pdf + static constexpr char qual_missing_code = -1; // 0xff from sam.c std::vector buffer; std::uint32_t name_len{}; @@ -40,9 +30,7 @@ class samrec { friend constexpr auto get_qual(const samrec &); friend constexpr auto get_qual_end(const samrec &); friend constexpr auto get_qual_size(const samrec &); - // clang-format on - // clang-format off samrec() = default; ~samrec() = default; samrec(const samrec &) = delete; @@ -52,37 +40,26 @@ class samrec { // clang-format on [[nodiscard]] auto - to_string() const { - // NOLINTBEGIN(cppcoreguidelines-pro-bounds-pointer-arithmetic) - return std::format( - "@{}\n{}\n+\n{}", - std::string(std::data(buffer), std::data(buffer) + name_len), - std::string(std::data(buffer) + name_len, - std::data(buffer) + name_len + seq_len), - std::string(std::data(buffer) + name_len + seq_len, - std::data(buffer) + name_len + seq_len + seq_len)); - // NOLINTEND(cppcoreguidelines-pro-bounds-pointer-arithmetic) - } + to_string() const -> std::string; operator bool() const { return name_len != 0; } [[nodiscard]] static auto - get_next(samrec::pos_t &cursor, const samrec::pos_t end_itr, - samrec &rec) -> bool; + get_next(pos_t &cursor, const pos_t end_itr, samrec &rec) -> bool; [[nodiscard]] static auto - find_end_pos(pos_t itr, const pos_t end) -> samrec::pos_t; + find_end_pos(pos_t itr, const pos_t end) -> pos_t; }; -// NOLINTBEGIN(cppcoreguidelines-pro-bounds-pointer-arithmetic) [[nodiscard]] inline constexpr auto get_name(const samrec &rec) { + // ADS: something is unhappy when this function returns an interator return std::data(rec.buffer); } [[nodiscard]] inline constexpr auto get_name_end(const samrec &rec) { - return get_name(rec) + rec.name_len; + return get_name(rec) + rec.name_len; // NOLINT } [[nodiscard]] inline constexpr auto @@ -92,7 +69,7 @@ get_seq(const samrec &rec) { [[nodiscard]] inline constexpr auto get_seq_end(const samrec &rec) { - return get_seq(rec) + rec.seq_len; + return get_seq(rec) + rec.seq_len; // NOLINT } [[nodiscard]] inline constexpr auto @@ -107,16 +84,13 @@ get_qual(const samrec &rec) { [[nodiscard]] inline constexpr auto get_qual_end(const samrec &rec) { - return get_qual(rec) + rec.seq_len; + return get_qual(rec) + rec.seq_len; // NOLINT } [[nodiscard]] inline constexpr auto get_qual_size(const samrec &rec) { return get_seq_size(rec); } -// NOLINTEND(cppcoreguidelines-pro-bounds-pointer-arithmetic) - -using sam_chunks_t = std::vector>; struct sam_task_t { samrec::pos_t beg{}; @@ -124,7 +98,7 @@ struct sam_task_t { }; [[nodiscard]] inline auto -get_next(auto &itr, const auto end, samrec &rec) -> bool { +get_next(samrec::pos_t &itr, const samrec::pos_t end, samrec &rec) -> bool { return samrec::get_next(itr, end, rec); } diff --git a/src/tile_processor.hpp b/src/tile_processor.hpp index 4f88e48..839144d 100644 --- a/src/tile_processor.hpp +++ b/src/tile_processor.hpp @@ -12,6 +12,7 @@ #include #include #include +#include #include #include #include @@ -103,7 +104,8 @@ class tile_processor { while (colon_count < tile_id_position && tile_itr != name_end) colon_count += (*tile_itr++ == ':'); std::uint32_t curr_tile_id{}; - const auto [_, ec] = std::from_chars(tile_itr, name_end, curr_tile_id); + const auto [_, ec] = std::from_chars( + std::to_address(tile_itr), std::to_address(name_end), curr_tile_id); if (ec != std::errc{}) throw std::system_error(std::make_error_code(ec), "failed to parse tile id");