diff --git a/CMakeLists.txt b/CMakeLists.txt index 03a5955..206812c 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -97,7 +97,7 @@ add_library(quality_score OBJECT src/quality_score.cpp) add_library(tile_processor OBJECT src/tile_processor.cpp) add_library(falco_grade OBJECT src/falco_grade.cpp) add_library(kmer_counter OBJECT src/kmer_counter.cpp) -add_library(contaminants OBJECT src/contaminants.cpp) +add_library(contaminant_set OBJECT src/contaminant_set.cpp) add_library(duplication_results OBJECT src/duplication_results.cpp) add_library(adapter_matcher OBJECT src/adapter_matcher.cpp) add_library(report OBJECT src/report.cpp) @@ -126,7 +126,7 @@ endif() add_executable(falco src/falco.cpp) target_link_libraries(falco PRIVATE tile_processor - contaminants + contaminant_set adapter_matcher adapter_set duplication_results diff --git a/src/bamrec.cpp b/src/bamrec.cpp index 287d169..b951e6c 100644 --- a/src/bamrec.cpp +++ b/src/bamrec.cpp @@ -13,15 +13,16 @@ [[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), // + 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); + auto qual_fixed = std::string(std::cbegin(qual_s), std::cend(qual_s)); + std::ranges::transform(qual_fixed, std::begin(qual_fixed), [&](const auto c) { + return c + quality_score_offset; + }); + 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)), qual_fixed); } @@ -33,33 +34,31 @@ bamrec::to_string() const -> std::string { #undef bam_seqi #endif -template -static inline constexpr OutputIt -assign_sequence_revcomp(BidirIt first, auto last, OutputIt d_first) { +template +static inline constexpr output_itr_t +assign_sequence_revcomp(bidir_itr_t first, auto last, output_itr_t d_first) { + static constexpr std::span seq_nt16_str = "=ACMGRSVTWYHKDBN"; 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"; +template +static inline constexpr output_itr_t +assign_sequence(bidir_itr_t first, auto last, output_itr_t d_first) { + static constexpr std::span 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; } diff --git a/src/contaminant_set.cpp b/src/contaminant_set.cpp new file mode 100644 index 0000000..64126d7 --- /dev/null +++ b/src/contaminant_set.cpp @@ -0,0 +1,274 @@ +// SPDX-License-Identifier: MIT; Copyright 2026 Andrew D Smith + +#include "contaminant_set.hpp" +#include "falco_utils.hpp" + +#include +#include +#include +#include +#include +#include +#include +#include +#include // for std::get +#include +#include + +// clang-format off +static constexpr auto default_contaminants = {std::pair + {R"(Illumina Single End Adapter 1)", R"(GATCGGAAGAGCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(Illumina Single End Adapter 2)", R"(CAAGCAGAAGACGGCATACGAGCTCTTCCGATCT)"}, + {R"(Illumina Single End PCR Primer 1)", R"(AATGATACGGCGACCACCGAGATCTACACTCTTTCCCTACACGACGCTCTTCCGATCT)"}, + {R"(Illumina Single End PCR Primer 2)", R"(CAAGCAGAAGACGGCATACGAGCTCTTCCGATCT)"}, + {R"(Illumina Single End Sequencing Primer)", R"(ACACTCTTTCCCTACACGACGCTCTTCCGATCT)"}, + {R"(Illumina Paired End Adapter 1)", R"(ACACTCTTTCCCTACACGACGCTCTTCCGATCT)"}, + {R"(Illumina Paired End Adapter 2)", R"(GATCGGAAGAGCGGTTCAGCAGGAATGCCGAG)"}, + {R"(Illumina Paried End PCR Primer 1)", R"(AATGATACGGCGACCACCGAGATCTACACTCTTTCCCTACACGACGCTCTTCCGATCT)"}, + {R"(Illumina Paired End PCR Primer 2)", R"(CAAGCAGAAGACGGCATACGAGATCGGTCTCGGCATTCCTGCTGAACCGCTCTTCCGATCT)"}, + {R"(Illumina Paried End Sequencing Primer 1)", R"(ACACTCTTTCCCTACACGACGCTCTTCCGATCT)"}, + {R"(Illumina Paired End Sequencing Primer 2)", R"(CGGTCTCGGCATTCCTGCTGAACCGCTCTTCCGATCT)"}, + {R"(Illumina DpnII expression Adapter 1)", R"(ACAGGTTCAGAGTTCTACAGTCCGAC)"}, + {R"(Illumina DpnII expression Adapter 2)", R"(CAAGCAGAAGACGGCATACGA)"}, + {R"(Illumina DpnII expression PCR Primer 1)", R"(CAAGCAGAAGACGGCATACGA)"}, + {R"(Illumina DpnII expression PCR Primer 2)", R"(AATGATACGGCGACCACCGACAGGTTCAGAGTTCTACAGTCCGA)"}, + {R"(Illumina DpnII expression Sequencing Primer)", R"(CGACAGGTTCAGAGTTCTACAGTCCGACGATC)"}, + {R"(Illumina NlaIII expression Adapter 1)", R"(ACAGGTTCAGAGTTCTACAGTCCGACATG)"}, + {R"(Illumina NlaIII expression Adapter 2)", R"(CAAGCAGAAGACGGCATACGA)"}, + {R"(Illumina NlaIII expression PCR Primer 1)", R"(CAAGCAGAAGACGGCATACGA)"}, + {R"(Illumina NlaIII expression PCR Primer 2)", R"(AATGATACGGCGACCACCGACAGGTTCAGAGTTCTACAGTCCGA)"}, + {R"(Illumina NlaIII expression Sequencing Primer)", R"(CCGACAGGTTCAGAGTTCTACAGTCCGACATG)"}, + {R"(Illumina Small RNA Adapter 1)", R"(GTTCAGAGTTCTACAGTCCGACGATC)"}, + {R"(Illumina Small RNA Adapter 2)", R"(TGGAATTCTCGGGTGCCAAGG)"}, + {R"(Illumina Small RNA RT Primer)", R"(CAAGCAGAAGACGGCATACGA)"}, + {R"(Illumina Small RNA PCR Primer 2)", R"(AATGATACGGCGACCACCGACAGGTTCAGAGTTCTACAGTCCGA)"}, + {R"(Illumina Small RNA Sequencing Primer)", R"(CGACAGGTTCAGAGTTCTACAGTCCGACGATC)"}, + {R"(Illumina Multiplexing Adapter 1)", R"(GATCGGAAGAGCACACGTCT)"}, + {R"(Illumina Multiplexing Adapter 2)", R"(ACACTCTTTCCCTACACGACGCTCTTCCGATCT)"}, + {R"(Illumina Multiplexing PCR Primer 1.01)", R"(AATGATACGGCGACCACCGAGATCTACACTCTTTCCCTACACGACGCTCTTCCGATCT)"}, + {R"(Illumina Multiplexing PCR Primer 2.01)", R"(GTGACTGGAGTTCAGACGTGTGCTCTTCCGATCT)"}, + {R"(Illumina Multiplexing Read1 Sequencing Primer)", R"(ACACTCTTTCCCTACACGACGCTCTTCCGATCT)"}, + {R"(Illumina Multiplexing Index Sequencing Primer)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCAC)"}, + {R"(Illumina Multiplexing Read2 Sequencing Primer)", R"(GTGACTGGAGTTCAGACGTGTGCTCTTCCGATCT)"}, + {R"(Illumina PCR Primer Index 1)", R"(CAAGCAGAAGACGGCATACGAGATCGTGATGTGACTGGAGTTC)"}, + {R"(Illumina PCR Primer Index 2)", R"(CAAGCAGAAGACGGCATACGAGATACATCGGTGACTGGAGTTC)"}, + {R"(Illumina PCR Primer Index 3)", R"(CAAGCAGAAGACGGCATACGAGATGCCTAAGTGACTGGAGTTC)"}, + {R"(Illumina PCR Primer Index 4)", R"(CAAGCAGAAGACGGCATACGAGATTGGTCAGTGACTGGAGTTC)"}, + {R"(Illumina PCR Primer Index 5)", R"(CAAGCAGAAGACGGCATACGAGATCACTGTGTGACTGGAGTTC)"}, + {R"(Illumina PCR Primer Index 6)", R"(CAAGCAGAAGACGGCATACGAGATATTGGCGTGACTGGAGTTC)"}, + {R"(Illumina PCR Primer Index 7)", R"(CAAGCAGAAGACGGCATACGAGATGATCTGGTGACTGGAGTTC)"}, + {R"(Illumina PCR Primer Index 8)", R"(CAAGCAGAAGACGGCATACGAGATTCAAGTGTGACTGGAGTTC)"}, + {R"(Illumina PCR Primer Index 9)", R"(CAAGCAGAAGACGGCATACGAGATCTGATCGTGACTGGAGTTC)"}, + {R"(Illumina PCR Primer Index 10)", R"(CAAGCAGAAGACGGCATACGAGATAAGCTAGTGACTGGAGTTC)"}, + {R"(Illumina PCR Primer Index 11)", R"(CAAGCAGAAGACGGCATACGAGATGTAGCCGTGACTGGAGTTC)"}, + {R"(Illumina PCR Primer Index 12)", R"(CAAGCAGAAGACGGCATACGAGATTACAAGGTGACTGGAGTTC)"}, + {R"(Illumina DpnII Gex Adapter 1)", R"(GATCGTCGGACTGTAGAACTCTGAAC)"}, + {R"(Illumina DpnII Gex Adapter 1.01)", R"(ACAGGTTCAGAGTTCTACAGTCCGAC)"}, + {R"(Illumina DpnII Gex Adapter 2)", R"(CAAGCAGAAGACGGCATACGA)"}, + {R"(Illumina DpnII Gex Adapter 2.01)", R"(TCGTATGCCGTCTTCTGCTTG)"}, + {R"(Illumina DpnII Gex PCR Primer 1)", R"(CAAGCAGAAGACGGCATACGA)"}, + {R"(Illumina DpnII Gex PCR Primer 2)", R"(AATGATACGGCGACCACCGACAGGTTCAGAGTTCTACAGTCCGA)"}, + {R"(Illumina DpnII Gex Sequencing Primer)", R"(CGACAGGTTCAGAGTTCTACAGTCCGACGATC)"}, + {R"(Illumina NlaIII Gex Adapter 1.01)", R"(TCGGACTGTAGAACTCTGAAC)"}, + {R"(Illumina NlaIII Gex Adapter 1.02)", R"(ACAGGTTCAGAGTTCTACAGTCCGACATG)"}, + {R"(Illumina NlaIII Gex Adapter 2.01)", R"(CAAGCAGAAGACGGCATACGA)"}, + {R"(Illumina NlaIII Gex Adapter 2.02)", R"(TCGTATGCCGTCTTCTGCTTG)"}, + {R"(Illumina NlaIII Gex PCR Primer 1)", R"(CAAGCAGAAGACGGCATACGA)"}, + {R"(Illumina NlaIII Gex PCR Primer 2)", R"(AATGATACGGCGACCACCGACAGGTTCAGAGTTCTACAGTCCGA)"}, + {R"(Illumina NlaIII Gex Sequencing Primer)", R"(CCGACAGGTTCAGAGTTCTACAGTCCGACATG)"}, + {R"(Illumina 5p RNA Adapter)", R"(GTTCAGAGTTCTACAGTCCGACGATC)"}, + {R"(Illumina RNA Adapter1)", R"(TGGAATTCTCGGGTGCCAAGG)"}, + {R"(Illumina Small RNA 3p Adapter 1)", R"(ATCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(Illumina Small RNA PCR Primer 1)", R"(CAAGCAGAAGACGGCATACGA)"}, + {R"(TruSeq Universal Adapter)", R"(AATGATACGGCGACCACCGAGATCTACACTCTTTCCCTACACGACGCTCTTCCGATCT)"}, + {R"(TruSeq Adapter, Index 1)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACATCACGATCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 2)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACCGATGTATCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 3)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACTTAGGCATCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 4)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACTGACCAATCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 5)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACACAGTGATCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 6)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACGCCAATATCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 7)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACCAGATCATCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 8)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACACTTGAATCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 9)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACGATCAGATCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 10)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACTAGCTTATCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 11)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACGGCTACATCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 12)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACCTTGTAATCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 13)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACAGTCAACTCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 14)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACAGTTCCGTCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 15)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACATGTCAGTCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 16)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACCCGTCCCTCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 18)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACGTCCGCATCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 19)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACGTGAAACTCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 20)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACGTGGCCTTCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 21)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACGTTTCGGTCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 22)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACCGTACGTTCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 23)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACCCACTCTTCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 25)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACACTGATATCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(TruSeq Adapter, Index 27)", R"(GATCGGAAGAGCACACGTCTGAACTCCAGTCACATTCCTTTCTCGTATGCCGTCTTCTGCTTG)"}, + {R"(Illumina RNA RT Primer)", R"(GCCTTGGCACCCGAGAATTCCA)"}, + {R"(Illumina RNA PCR Primer)", R"(AATGATACGGCGACCACCGAGATCTACACGTTCAGAGTTCTACAGTCCGA)"}, + {R"(RNA PCR Primer, Index 1)", R"(CAAGCAGAAGACGGCATACGAGATCGTGATGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 2)", R"(CAAGCAGAAGACGGCATACGAGATACATCGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 3)", R"(CAAGCAGAAGACGGCATACGAGATGCCTAAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 4)", R"(CAAGCAGAAGACGGCATACGAGATTGGTCAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 5)", R"(CAAGCAGAAGACGGCATACGAGATCACTGTGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 6)", R"(CAAGCAGAAGACGGCATACGAGATATTGGCGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 7)", R"(CAAGCAGAAGACGGCATACGAGATGATCTGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 8)", R"(CAAGCAGAAGACGGCATACGAGATTCAAGTGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 9)", R"(CAAGCAGAAGACGGCATACGAGATCTGATCGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 10)", R"(CAAGCAGAAGACGGCATACGAGATAAGCTAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 11)", R"(CAAGCAGAAGACGGCATACGAGATGTAGCCGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 12)", R"(CAAGCAGAAGACGGCATACGAGATTACAAGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 13)", R"(CAAGCAGAAGACGGCATACGAGATTTGACTGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 14)", R"(CAAGCAGAAGACGGCATACGAGATGGAACTGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 15)", R"(CAAGCAGAAGACGGCATACGAGATTGACATGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 16)", R"(CAAGCAGAAGACGGCATACGAGATGGACGGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 17)", R"(CAAGCAGAAGACGGCATACGAGATCTCTACGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 18)", R"(CAAGCAGAAGACGGCATACGAGATGCGGACGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 19)", R"(CAAGCAGAAGACGGCATACGAGATTTTCACGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 20)", R"(CAAGCAGAAGACGGCATACGAGATGGCCACGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 21)", R"(CAAGCAGAAGACGGCATACGAGATCGAAACGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 22)", R"(CAAGCAGAAGACGGCATACGAGATCGTACGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 23)", R"(CAAGCAGAAGACGGCATACGAGATCCACTCGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 24)", R"(CAAGCAGAAGACGGCATACGAGATGCTACCGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 25)", R"(CAAGCAGAAGACGGCATACGAGATATCAGTGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 26)", R"(CAAGCAGAAGACGGCATACGAGATGCTCATGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 27)", R"(CAAGCAGAAGACGGCATACGAGATAGGAATGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 28)", R"(CAAGCAGAAGACGGCATACGAGATCTTTTGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 29)", R"(CAAGCAGAAGACGGCATACGAGATTAGTTGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 30)", R"(CAAGCAGAAGACGGCATACGAGATCCGGTGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 31)", R"(CAAGCAGAAGACGGCATACGAGATATCGTGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 32)", R"(CAAGCAGAAGACGGCATACGAGATTGAGTGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 33)", R"(CAAGCAGAAGACGGCATACGAGATCGCCTGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 34)", R"(CAAGCAGAAGACGGCATACGAGATGCCATGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 35)", R"(CAAGCAGAAGACGGCATACGAGATAAAATGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 36)", R"(CAAGCAGAAGACGGCATACGAGATTGTTGGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 37)", R"(CAAGCAGAAGACGGCATACGAGATATTCCGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 38)", R"(CAAGCAGAAGACGGCATACGAGATAGCTAGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 39)", R"(CAAGCAGAAGACGGCATACGAGATGTATAGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 40)", R"(CAAGCAGAAGACGGCATACGAGATTCTGAGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 41)", R"(CAAGCAGAAGACGGCATACGAGATGTCGTCGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 42)", R"(CAAGCAGAAGACGGCATACGAGATCGATTAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 43)", R"(CAAGCAGAAGACGGCATACGAGATGCTGTAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 44)", R"(CAAGCAGAAGACGGCATACGAGATATTATAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 45)", R"(CAAGCAGAAGACGGCATACGAGATGAATGAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 46)", R"(CAAGCAGAAGACGGCATACGAGATTCGGGAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 47)", R"(CAAGCAGAAGACGGCATACGAGATCTTCGAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(RNA PCR Primer, Index 48)", R"(CAAGCAGAAGACGGCATACGAGATTGCCGAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA)"}, + {R"(ABI Dynabead EcoP Oligo)", R"(CTGATCTAGAGGTACCGGATCCCAGCAGT)"}, + {R"(ABI Solid3 Adapter A)", R"(CTGCCCCGGGTTCCTCATTCTCTCAGCAGCATG)"}, + {R"(ABI Solid3 Adapter B)", R"(CCACTACGCCTCCGCTTTCCTCTCTATGGGCAGTCGGTGAT)"}, + {R"(ABI Solid3 5' AMP Primer)", R"(CCACTACGCCTCCGCTTTCCTCTCTATG)"}, + {R"(ABI Solid3 3' AMP Primer)", R"(CTGCCCCGGGTTCCTCATTCT)"}, + {R"(ABI Solid3 EF1 alpha Sense Primer)", R"(CATGTGTGTTGAGAGCTTC)"}, + {R"(ABI Solid3 EF1 alpha Antisense Primer)", R"(GAAAACCAAAGTGGTCCAC)"}, + {R"(ABI Solid3 GAPDH Forward Primer)", R"(TTAGCACCCCTGGCCAAGG)"}, + {R"(ABI Solid3 GAPDH Reverse Primer)", R"(CTTACTCCTTGGAGGCCATG)"}, + {R"(Clontech Universal Primer Mix Short)", R"(CTAATACGACTCACTATAGGGC)"}, + {R"(Clontech Universal Primer Mix Long)", R"(CTAATACGACTCACTATAGGGCAAGCAGTGGTATCAACGCAGAGT)"}, + {R"(Clontech SMARTer II A Oligonucleotide)", R"(AAGCAGTGGTATCAACGCAGAGTAC)"}, + {R"(Clontech SMART CDS Primer II A)", R"(AAGCAGTGGTATCAACGCAGAGTACT)"}, + {R"(Clontech_Universal_Primer_Mix_Short)", R"(CTAATACGACTCACTATAGGGC)"}, + {R"(Clontech_Universal_Primer_Mix_Long)", R"(CTAATACGACTCACTATAGGGCAAGCAGTGGTATCAACGCAGAGT)"}, + {R"(Clontech_SMARTer_II_A_Oligonucleotide)", R"(AAGCAGTGGTATCAACGCAGAGTAC)"}, + {R"(Clontech_SMART_CDS_Primer_II_A)", R"(AAGCAGTGGTATCAACGCAGAGTACT)"}, + {R"(Clontech_SMART_CDS_Primer_II_A)", R"(ACGTACTCTGCGTTGATACCACTGCTTCCGCGGACAGGCGTGTAGATCTCGGTGGTCGC)"}, + {R"(Clontech_SMART_CDS_Primer_II_A)", R"(GAGTACGTACTCTGCGTTGATACCACTGCTTCCGCGGACAGGCGTGTAGATCTCGGTGGT)"}, +}; +// clang-format on + +[[nodiscard]] static auto +load_contaminants(const std::string &filename) + -> std::vector> { + // ADS: (todo) handle carriage returns and other control chars + std::ifstream in(filename); + if (!in) + throw std::runtime_error("failed to open contaminants file: " + filename); + std::vector> contaminants; + std::string line_data; + while (std::getline(in, line_data)) { + std::string_view line = line_data; + const auto to_keep_prefix = line.find_first_not_of(" \t"); + if (to_keep_prefix == std::string_view::npos) + continue; + line.remove_prefix(std::min(to_keep_prefix, std::size(line))); + if (line[0] == '#') + continue; + const auto to_keep_suffix = line.find_last_not_of(" \t"); + if (to_keep_suffix == std::string_view::npos) + continue; + line.remove_suffix(std::size(line) - to_keep_suffix - 1); + std::string cleaned_line; + for (auto itr = std::cbegin(line); itr != std::cend(line); + itr = std::next(itr)) + if (!std::isblank(*itr) || + (std::next(itr) != std::cend(line) && *itr != *std::next(itr))) + cleaned_line += *itr; + const auto tab_pos = cleaned_line.find('\t'); + if (tab_pos == std::string::npos || + tab_pos != cleaned_line.find_last_of('\t')) + throw std::runtime_error("malformed line: " + line_data); + const auto is_print = [](const auto c) { return std::isprint(c); }; + const auto name = cleaned_line.substr(0, tab_pos); + const auto seq = cleaned_line.substr(tab_pos + 1); + if (!std::ranges::all_of(name, is_print) || + !std::ranges::all_of(seq, is_print)) + throw std::runtime_error("malformed line: " + line_data); + contaminants.emplace_back(name, seq); + } + return contaminants; +} + +[[nodiscard]] auto +contaminant_set::get_name(const std::int64_t contam_idx) + -> const std::string & { + using std::string_literals::operator""s; + static constexpr auto no_hit_label = "No Hit"s; + const auto &contaminants = instance().contaminants; + if (contam_idx < 0 || contam_idx >= std::ssize(contaminants)) + return no_hit_label; + return contaminants[contam_idx].first; +} + +/// get the longest substring of left that is a prefix of right +[[nodiscard]] static inline auto +get_overlap(const auto &left, const auto &right) { + const auto left_beg = std::cbegin(left); + const auto left_end = std::cend(left); + auto best_n_matches = 0L; + for (auto left_itr = left_beg; left_itr != left_end; ++left_itr) { + const auto [mm_left, _] = + std::ranges::mismatch(std::ranges::subrange(left_itr, left_end), right); + const auto n_matches = std::distance(left_itr, mm_left); + best_n_matches = std::max(best_n_matches, n_matches); + } + return best_n_matches; +} + +[[nodiscard]] auto +contaminant_set::match(const std::string &query) -> std::int64_t { + const auto &contaminants = instance().contaminants; + auto best_idx = 0L; + auto best_match = 0L; + auto best_match_len = 0L; + for (const auto &[idx, seq] : + falco::views::enumerate(std::views::elements<1>(contaminants))) { + const auto n_match = + std::max(get_overlap(query, seq), get_overlap(seq, query)); + if (n_match > best_match) { + best_idx = idx; + best_match = n_match; + best_match_len = std::ssize(seq); + } + } + const auto match_cutoff = std::min(best_match_len, std::ssize(query)) / 2; + // If any sequence is a match, return the best one + return best_match < match_cutoff ? -1 : best_idx; +} + +contaminant_set::contaminant_set(const std::string &filename) { + if (filename.empty()) + std::ranges::copy(default_contaminants, std::back_inserter(contaminants)); + else + contaminants = load_contaminants(filename); +} diff --git a/src/contaminant_set.hpp b/src/contaminant_set.hpp new file mode 100644 index 0000000..bd38475 --- /dev/null +++ b/src/contaminant_set.hpp @@ -0,0 +1,45 @@ +// SPDX-License-Identifier: MIT; Copyright 2026 Andrew D Smith + +#ifndef SRC_CONTAMINANT_SET_HPP_ +#define SRC_CONTAMINANT_SET_HPP_ + +#include +#include +#include // for std::pair +#include +#include +#include + +struct contaminant_set { + static auto + instance(const std::string &filename = std::string{}) + -> const contaminant_set & { + static const contaminant_set s(filename); + return s; + } + + [[nodiscard]] static auto + n_contaminants() -> std::uint64_t { + return std::size(instance().contaminants); + } + + [[nodiscard]] static auto + get_name(std::int64_t idx) -> const std::string &; + + [[nodiscard]] static auto + match(const std::string &query) -> std::int64_t; + + // clang-format off + contaminant_set(const contaminant_set &) = delete; + contaminant_set(contaminant_set &&) = delete; + auto operator=(const contaminant_set &) = delete; + auto operator=(contaminant_set &&) = delete; + ~contaminant_set() = default; + // clang-format on + +private: + std::vector> contaminants; + explicit contaminant_set(const std::string &filename); +}; // contaminant_set + +#endif // SRC_CONTAMINANT_SET_HPP_ diff --git a/src/contaminants.cpp b/src/contaminants.cpp deleted file mode 100644 index d6b64ef..0000000 --- a/src/contaminants.cpp +++ /dev/null @@ -1,262 +0,0 @@ -// SPDX-License-Identifier: MIT; Copyright 2026 Andrew D Smith - -#include "contaminants.hpp" -#include "falco_utils.hpp" - -#include -#include -#include -#include -#include -#include -#include -#include -#include // for std::get -#include -#include - -// clang-format off -std::vector> contaminants = { // NOLINT(cert-err58-cpp,cppcoreguidelines-avoid-non-const-global-variables) - {"Illumina Single End Adapter 1", "GATCGGAAGAGCTCGTATGCCGTCTTCTGCTTG"}, - {"Illumina Single End Adapter 2", "CAAGCAGAAGACGGCATACGAGCTCTTCCGATCT"}, - {"Illumina Single End PCR Primer 1", "AATGATACGGCGACCACCGAGATCTACACTCTTTCCCTACACGACGCTCTTCCGATCT"}, - {"Illumina Single End PCR Primer 2", "CAAGCAGAAGACGGCATACGAGCTCTTCCGATCT"}, - {"Illumina Single End Sequencing Primer", "ACACTCTTTCCCTACACGACGCTCTTCCGATCT"}, - {"Illumina Paired End Adapter 1", "ACACTCTTTCCCTACACGACGCTCTTCCGATCT"}, - {"Illumina Paired End Adapter 2", "GATCGGAAGAGCGGTTCAGCAGGAATGCCGAG"}, - {"Illumina Paried End PCR Primer 1", "AATGATACGGCGACCACCGAGATCTACACTCTTTCCCTACACGACGCTCTTCCGATCT"}, - {"Illumina Paired End PCR Primer 2", "CAAGCAGAAGACGGCATACGAGATCGGTCTCGGCATTCCTGCTGAACCGCTCTTCCGATCT"}, - {"Illumina Paried End Sequencing Primer 1", "ACACTCTTTCCCTACACGACGCTCTTCCGATCT"}, - {"Illumina Paired End Sequencing Primer 2", "CGGTCTCGGCATTCCTGCTGAACCGCTCTTCCGATCT"}, - {"Illumina DpnII expression Adapter 1", "ACAGGTTCAGAGTTCTACAGTCCGAC"}, - {"Illumina DpnII expression Adapter 2", "CAAGCAGAAGACGGCATACGA"}, - {"Illumina DpnII expression PCR Primer 1", "CAAGCAGAAGACGGCATACGA"}, - {"Illumina DpnII expression PCR Primer 2", "AATGATACGGCGACCACCGACAGGTTCAGAGTTCTACAGTCCGA"}, - {"Illumina DpnII expression Sequencing Primer", "CGACAGGTTCAGAGTTCTACAGTCCGACGATC"}, - {"Illumina NlaIII expression Adapter 1", "ACAGGTTCAGAGTTCTACAGTCCGACATG"}, - {"Illumina NlaIII expression Adapter 2", "CAAGCAGAAGACGGCATACGA"}, - {"Illumina NlaIII expression PCR Primer 1", "CAAGCAGAAGACGGCATACGA"}, - {"Illumina NlaIII expression PCR Primer 2", "AATGATACGGCGACCACCGACAGGTTCAGAGTTCTACAGTCCGA"}, - {"Illumina NlaIII expression Sequencing Primer", "CCGACAGGTTCAGAGTTCTACAGTCCGACATG"}, - {"Illumina Small RNA Adapter 1", "GTTCAGAGTTCTACAGTCCGACGATC"}, - {"Illumina Small RNA Adapter 2", "TGGAATTCTCGGGTGCCAAGG"}, - {"Illumina Small RNA RT Primer", "CAAGCAGAAGACGGCATACGA"}, - {"Illumina Small RNA PCR Primer 2", "AATGATACGGCGACCACCGACAGGTTCAGAGTTCTACAGTCCGA"}, - {"Illumina Small RNA Sequencing Primer", "CGACAGGTTCAGAGTTCTACAGTCCGACGATC"}, - {"Illumina Multiplexing Adapter 1", "GATCGGAAGAGCACACGTCT"}, - {"Illumina Multiplexing Adapter 2", "ACACTCTTTCCCTACACGACGCTCTTCCGATCT"}, - {"Illumina Multiplexing PCR Primer 1.01", "AATGATACGGCGACCACCGAGATCTACACTCTTTCCCTACACGACGCTCTTCCGATCT"}, - {"Illumina Multiplexing PCR Primer 2.01", "GTGACTGGAGTTCAGACGTGTGCTCTTCCGATCT"}, - {"Illumina Multiplexing Read1 Sequencing Primer", "ACACTCTTTCCCTACACGACGCTCTTCCGATCT"}, - {"Illumina Multiplexing Index Sequencing Primer", "GATCGGAAGAGCACACGTCTGAACTCCAGTCAC"}, - {"Illumina Multiplexing Read2 Sequencing Primer", "GTGACTGGAGTTCAGACGTGTGCTCTTCCGATCT"}, - {"Illumina PCR Primer Index 1", "CAAGCAGAAGACGGCATACGAGATCGTGATGTGACTGGAGTTC"}, - {"Illumina PCR Primer Index 2", "CAAGCAGAAGACGGCATACGAGATACATCGGTGACTGGAGTTC"}, - {"Illumina PCR Primer Index 3", "CAAGCAGAAGACGGCATACGAGATGCCTAAGTGACTGGAGTTC"}, - {"Illumina PCR Primer Index 4", "CAAGCAGAAGACGGCATACGAGATTGGTCAGTGACTGGAGTTC"}, - {"Illumina PCR Primer Index 5", "CAAGCAGAAGACGGCATACGAGATCACTGTGTGACTGGAGTTC"}, - {"Illumina PCR Primer Index 6", "CAAGCAGAAGACGGCATACGAGATATTGGCGTGACTGGAGTTC"}, - {"Illumina PCR Primer Index 7", "CAAGCAGAAGACGGCATACGAGATGATCTGGTGACTGGAGTTC"}, - {"Illumina PCR Primer Index 8", "CAAGCAGAAGACGGCATACGAGATTCAAGTGTGACTGGAGTTC"}, - {"Illumina PCR Primer Index 9", "CAAGCAGAAGACGGCATACGAGATCTGATCGTGACTGGAGTTC"}, - {"Illumina PCR Primer Index 10", "CAAGCAGAAGACGGCATACGAGATAAGCTAGTGACTGGAGTTC"}, - {"Illumina PCR Primer Index 11", "CAAGCAGAAGACGGCATACGAGATGTAGCCGTGACTGGAGTTC"}, - {"Illumina PCR Primer Index 12", "CAAGCAGAAGACGGCATACGAGATTACAAGGTGACTGGAGTTC"}, - {"Illumina DpnII Gex Adapter 1", "GATCGTCGGACTGTAGAACTCTGAAC"}, - {"Illumina DpnII Gex Adapter 1.01", "ACAGGTTCAGAGTTCTACAGTCCGAC"}, - {"Illumina DpnII Gex Adapter 2", "CAAGCAGAAGACGGCATACGA"}, - {"Illumina DpnII Gex Adapter 2.01", "TCGTATGCCGTCTTCTGCTTG"}, - {"Illumina DpnII Gex PCR Primer 1", "CAAGCAGAAGACGGCATACGA"}, - {"Illumina DpnII Gex PCR Primer 2", "AATGATACGGCGACCACCGACAGGTTCAGAGTTCTACAGTCCGA"}, - {"Illumina DpnII Gex Sequencing Primer", "CGACAGGTTCAGAGTTCTACAGTCCGACGATC"}, - {"Illumina NlaIII Gex Adapter 1.01", "TCGGACTGTAGAACTCTGAAC"}, - {"Illumina NlaIII Gex Adapter 1.02", "ACAGGTTCAGAGTTCTACAGTCCGACATG"}, - {"Illumina NlaIII Gex Adapter 2.01", "CAAGCAGAAGACGGCATACGA"}, - {"Illumina NlaIII Gex Adapter 2.02", "TCGTATGCCGTCTTCTGCTTG"}, - {"Illumina NlaIII Gex PCR Primer 1", "CAAGCAGAAGACGGCATACGA"}, - {"Illumina NlaIII Gex PCR Primer 2", "AATGATACGGCGACCACCGACAGGTTCAGAGTTCTACAGTCCGA"}, - {"Illumina NlaIII Gex Sequencing Primer", "CCGACAGGTTCAGAGTTCTACAGTCCGACATG"}, - {"Illumina 5p RNA Adapter", "GTTCAGAGTTCTACAGTCCGACGATC"}, - {"Illumina RNA Adapter1", "TGGAATTCTCGGGTGCCAAGG"}, - {"Illumina Small RNA 3p Adapter 1", "ATCTCGTATGCCGTCTTCTGCTTG"}, - {"Illumina Small RNA PCR Primer 1", "CAAGCAGAAGACGGCATACGA"}, - {"TruSeq Universal Adapter", "AATGATACGGCGACCACCGAGATCTACACTCTTTCCCTACACGACGCTCTTCCGATCT"}, - {"TruSeq Adapter, Index 1", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACATCACGATCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 2", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACCGATGTATCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 3", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACTTAGGCATCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 4", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACTGACCAATCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 5", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACACAGTGATCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 6", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACGCCAATATCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 7", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACCAGATCATCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 8", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACACTTGAATCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 9", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACGATCAGATCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 10", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACTAGCTTATCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 11", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACGGCTACATCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 12", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACCTTGTAATCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 13", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACAGTCAACTCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 14", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACAGTTCCGTCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 15", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACATGTCAGTCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 16", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACCCGTCCCTCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 18", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACGTCCGCATCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 19", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACGTGAAACTCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 20", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACGTGGCCTTCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 21", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACGTTTCGGTCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 22", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACCGTACGTTCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 23", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACCCACTCTTCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 25", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACACTGATATCTCGTATGCCGTCTTCTGCTTG"}, - {"TruSeq Adapter, Index 27", "GATCGGAAGAGCACACGTCTGAACTCCAGTCACATTCCTTTCTCGTATGCCGTCTTCTGCTTG"}, - {"Illumina RNA RT Primer", "GCCTTGGCACCCGAGAATTCCA"}, - {"Illumina RNA PCR Primer", "AATGATACGGCGACCACCGAGATCTACACGTTCAGAGTTCTACAGTCCGA"}, - {"RNA PCR Primer, Index 1", "CAAGCAGAAGACGGCATACGAGATCGTGATGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 2", "CAAGCAGAAGACGGCATACGAGATACATCGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 3", "CAAGCAGAAGACGGCATACGAGATGCCTAAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 4", "CAAGCAGAAGACGGCATACGAGATTGGTCAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 5", "CAAGCAGAAGACGGCATACGAGATCACTGTGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 6", "CAAGCAGAAGACGGCATACGAGATATTGGCGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 7", "CAAGCAGAAGACGGCATACGAGATGATCTGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 8", "CAAGCAGAAGACGGCATACGAGATTCAAGTGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 9", "CAAGCAGAAGACGGCATACGAGATCTGATCGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 10", "CAAGCAGAAGACGGCATACGAGATAAGCTAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 11", "CAAGCAGAAGACGGCATACGAGATGTAGCCGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 12", "CAAGCAGAAGACGGCATACGAGATTACAAGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 13", "CAAGCAGAAGACGGCATACGAGATTTGACTGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 14", "CAAGCAGAAGACGGCATACGAGATGGAACTGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 15", "CAAGCAGAAGACGGCATACGAGATTGACATGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 16", "CAAGCAGAAGACGGCATACGAGATGGACGGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 17", "CAAGCAGAAGACGGCATACGAGATCTCTACGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 18", "CAAGCAGAAGACGGCATACGAGATGCGGACGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 19", "CAAGCAGAAGACGGCATACGAGATTTTCACGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 20", "CAAGCAGAAGACGGCATACGAGATGGCCACGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 21", "CAAGCAGAAGACGGCATACGAGATCGAAACGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 22", "CAAGCAGAAGACGGCATACGAGATCGTACGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 23", "CAAGCAGAAGACGGCATACGAGATCCACTCGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 24", "CAAGCAGAAGACGGCATACGAGATGCTACCGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 25", "CAAGCAGAAGACGGCATACGAGATATCAGTGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 26", "CAAGCAGAAGACGGCATACGAGATGCTCATGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 27", "CAAGCAGAAGACGGCATACGAGATAGGAATGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 28", "CAAGCAGAAGACGGCATACGAGATCTTTTGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 29", "CAAGCAGAAGACGGCATACGAGATTAGTTGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 30", "CAAGCAGAAGACGGCATACGAGATCCGGTGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 31", "CAAGCAGAAGACGGCATACGAGATATCGTGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 32", "CAAGCAGAAGACGGCATACGAGATTGAGTGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 33", "CAAGCAGAAGACGGCATACGAGATCGCCTGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 34", "CAAGCAGAAGACGGCATACGAGATGCCATGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 35", "CAAGCAGAAGACGGCATACGAGATAAAATGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 36", "CAAGCAGAAGACGGCATACGAGATTGTTGGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 37", "CAAGCAGAAGACGGCATACGAGATATTCCGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 38", "CAAGCAGAAGACGGCATACGAGATAGCTAGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 39", "CAAGCAGAAGACGGCATACGAGATGTATAGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 40", "CAAGCAGAAGACGGCATACGAGATTCTGAGGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 41", "CAAGCAGAAGACGGCATACGAGATGTCGTCGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 42", "CAAGCAGAAGACGGCATACGAGATCGATTAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 43", "CAAGCAGAAGACGGCATACGAGATGCTGTAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 44", "CAAGCAGAAGACGGCATACGAGATATTATAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 45", "CAAGCAGAAGACGGCATACGAGATGAATGAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 46", "CAAGCAGAAGACGGCATACGAGATTCGGGAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 47", "CAAGCAGAAGACGGCATACGAGATCTTCGAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"RNA PCR Primer, Index 48", "CAAGCAGAAGACGGCATACGAGATTGCCGAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA"}, - {"ABI Dynabead EcoP Oligo", "CTGATCTAGAGGTACCGGATCCCAGCAGT"}, - {"ABI Solid3 Adapter A", "CTGCCCCGGGTTCCTCATTCTCTCAGCAGCATG"}, - {"ABI Solid3 Adapter B", "CCACTACGCCTCCGCTTTCCTCTCTATGGGCAGTCGGTGAT"}, - {"ABI Solid3 5' AMP Primer", "CCACTACGCCTCCGCTTTCCTCTCTATG"}, - {"ABI Solid3 3' AMP Primer", "CTGCCCCGGGTTCCTCATTCT"}, - {"ABI Solid3 EF1 alpha Sense Primer", "CATGTGTGTTGAGAGCTTC"}, - {"ABI Solid3 EF1 alpha Antisense Primer", "GAAAACCAAAGTGGTCCAC"}, - {"ABI Solid3 GAPDH Forward Primer", "TTAGCACCCCTGGCCAAGG"}, - {"ABI Solid3 GAPDH Reverse Primer", "CTTACTCCTTGGAGGCCATG"}, - {"Clontech Universal Primer Mix Short", "CTAATACGACTCACTATAGGGC"}, - {"Clontech Universal Primer Mix Long", "CTAATACGACTCACTATAGGGCAAGCAGTGGTATCAACGCAGAGT"}, - {"Clontech SMARTer II A Oligonucleotide", "AAGCAGTGGTATCAACGCAGAGTAC"}, - {"Clontech SMART CDS Primer II A", "AAGCAGTGGTATCAACGCAGAGTACT"}, - {"Clontech_Universal_Primer_Mix_Short", "CTAATACGACTCACTATAGGGC"}, - {"Clontech_Universal_Primer_Mix_Long", "CTAATACGACTCACTATAGGGCAAGCAGTGGTATCAACGCAGAGT"}, - {"Clontech_SMARTer_II_A_Oligonucleotide", "AAGCAGTGGTATCAACGCAGAGTAC"}, - {"Clontech_SMART_CDS_Primer_II_A", "AAGCAGTGGTATCAACGCAGAGTACT"}, - {"Clontech_SMART_CDS_Primer_II_A", "ACGTACTCTGCGTTGATACCACTGCTTCCGCGGACAGGCGTGTAGATCTCGGTGGTCGC"}, - {"Clontech_SMART_CDS_Primer_II_A", "GAGTACGTACTCTGCGTTGATACCACTGCTTCCGCGGACAGGCGTGTAGATCTCGGTGGT"}, -}; -// clang-format on - -auto -load_contaminants(const std::string &filename) -> void { - // ADS: (todo) handle carriage returns and other control chars - std::ifstream in(filename); - if (!in) - throw std::runtime_error("failed to open contaminants file: " + filename); - contaminants.clear(); - std::string line_data; - while (std::getline(in, line_data)) { - std::string_view line = line_data; - const auto to_keep_prefix = line.find_first_not_of(" \t"); - if (to_keep_prefix == std::string_view::npos) - continue; - line.remove_prefix(std::min(to_keep_prefix, std::size(line))); - if (line[0] == '#') - continue; - const auto to_keep_suffix = line.find_last_not_of(" \t"); - if (to_keep_suffix == std::string_view::npos) - continue; - line.remove_suffix(std::size(line) - to_keep_suffix - 1); - std::string cleaned_line; - for (auto itr = std::cbegin(line); itr != std::cend(line); - itr = std::next(itr)) - if (!std::isblank(*itr) || - (std::next(itr) != std::cend(line) && *itr != *std::next(itr))) - cleaned_line += *itr; - const auto tab_pos = cleaned_line.find('\t'); - if (tab_pos == std::string::npos || - tab_pos != cleaned_line.find_last_of('\t')) - throw std::runtime_error("malformed line: " + line_data); - const auto is_print = [](const auto c) { return std::isprint(c); }; - const auto name = cleaned_line.substr(0, tab_pos); - const auto seq = cleaned_line.substr(tab_pos + 1); - if (!std::ranges::all_of(name, is_print) || - !std::ranges::all_of(seq, is_print)) - throw std::runtime_error("malformed line: " + line_data); - contaminants.emplace_back(name, seq); - } -} - -[[nodiscard]] auto -get_contam_name(const std::int64_t contam_idx) -> const std::string & { - using std::string_literals::operator""s; - static constexpr auto no_hit_label = "No Hit"s; - if (contam_idx < 0 || contam_idx >= std::ssize(contaminants)) - return no_hit_label; - return contaminants[contam_idx].first; -} - -// get the longest substring of left that is a prefix of right -[[nodiscard]] static inline auto -get_overlap(const auto &left, const auto &right) { - const auto left_beg = std::cbegin(left); - const auto left_end = std::cend(left); - auto best_n_matches = 0L; - for (auto left_itr = left_beg; left_itr != left_end; ++left_itr) { - const auto [mm_left, _] = - std::ranges::mismatch(std::ranges::subrange(left_itr, left_end), right); - const auto n_matches = std::distance(left_itr, mm_left); - best_n_matches = std::max(best_n_matches, n_matches); - } - return best_n_matches; -} - -[[nodiscard]] auto -match_contaminant(const std::string &query) -> std::int64_t { - auto best_idx = 0L; - auto best_match = 0L; - auto best_match_len = 0L; - for (const auto &[idx, seq] : - falco::views::enumerate(std::views::elements<1>(contaminants))) { - const auto n_match = - std::max(get_overlap(query, seq), get_overlap(seq, query)); - if (n_match > best_match) { - best_idx = idx; - best_match = n_match; - best_match_len = std::ssize(seq); - } - } - const auto match_cutoff = std::min(best_match_len, std::ssize(query)) / 2; - // If any sequence is a match, return the best one - return best_match < match_cutoff ? -1 : best_idx; -} diff --git a/src/contaminants.hpp b/src/contaminants.hpp deleted file mode 100644 index fba5a64..0000000 --- a/src/contaminants.hpp +++ /dev/null @@ -1,24 +0,0 @@ -// SPDX-License-Identifier: MIT; Copyright 2026 Andrew D Smith - -#ifndef SRC_CONTAMINANTS_HPP_ -#define SRC_CONTAMINANTS_HPP_ - -#include -#include // for std::pair -#include -#include // IWYU pragma: keep -#include - -// NOLINTNEXTLINE(cppcoreguidelines-avoid-non-const-global-variables) -extern std::vector> contaminants; - -[[nodiscard]] auto -get_contam_name(std::int64_t contam_idx) -> const std::string &; - -auto -load_contaminants(const std::string &filename) -> void; - -[[nodiscard]] auto -match_contaminant(const std::string &query) -> std::int64_t; - -#endif // SRC_CONTAMINANTS_HPP_ diff --git a/src/duplication_results.cpp b/src/duplication_results.cpp index 05e79d8..a1f6a03 100644 --- a/src/duplication_results.cpp +++ b/src/duplication_results.cpp @@ -2,7 +2,7 @@ #include "duplication_results.hpp" -#include "contaminants.hpp" +#include "contaminant_set.hpp" #include "falco_grade.hpp" #include "falco_utils.hpp" #include "falco_word.hpp" @@ -122,13 +122,13 @@ duplication_results::get_overrepresented(const std::uint64_t n_reads) const ret.emplace_back( seq, n_obs, pct(as_frac(n_obs, std::max(static_cast(1), n_reads))), - match_contaminant(seq.string())); + contaminant_set::match(seq.string())); return ret; } auto -duplication_results::initialize(const run_mode &mode, const file_info &info) - -> void { +duplication_results::initialize(const run_mode &mode, + const file_info &info) -> void { read_skip = info.n_reads_est < max_n_reads_total ? 0 @@ -177,7 +177,7 @@ overrepresented_report(const std::vector &overrep, r += header; for (const auto &[seq, n_obs, pct_val, contam_id] : overrep) r += std::format("{}\t{}\t{:.3g}\t{}\n", seq, n_obs, pct_val, - get_contam_name(contam_id)); + contaminant_set::get_name(contam_id)); } return r + end_module_tag; } @@ -285,8 +285,8 @@ get_grade_duplication(const dup_summary_t &summary) -> std::string { } [[nodiscard]] auto -duplication_report(const dup_summary_t &summary, const file_grades &grades) - -> std::string { +duplication_report(const dup_summary_t &summary, + const file_grades &grades) -> std::string { static constexpr auto label = "duplication"; static constexpr auto start_tag = ">>Sequence Duplication Levels\t{}\n" "#Total Deduplicated Percentage\t{:.6f}\n"; @@ -342,15 +342,15 @@ overrepresented_html(const std::vector &overrep, "No overrepresented sequences"); const auto rows = std::views::transform(overrep, [&](const auto &o) { return fmt::format(html_table_row_fmt, o.w.string(), o.n_obs, o.pct_val, - get_contam_name(o.contam_id)); + contaminant_set::get_name(o.contam_id)); }); return fmt::format(html_module_fmt, grade, label, title, grade, fmt::format(html_table, fmt::join(rows, "\n"))); } [[nodiscard]] auto -duplication_html(const dup_summary_t &summary, const file_grades &grades) - -> std::string { +duplication_html(const dup_summary_t &summary, + const file_grades &grades) -> std::string { static constexpr auto label = "duplication"; static constexpr auto plot_format = R"(