diff --git a/src/config/option_collation.cpp b/src/config/option_collation.cpp index 93d12efb2..0d819363d 100644 --- a/src/config/option_collation.cpp +++ b/src/config/option_collation.cpp @@ -2061,6 +2061,11 @@ boost::optional split_bamout_request(const OptionMap& options) return boost::none; } +bool full_bamouts_requested(const OptionMap& options) +{ + return options.at("full-bamout").as(); +} + unsigned max_open_read_files(const OptionMap& options) { return 2 * std::min(as_unsigned("max-open-read-files", options), count_read_paths(options)); diff --git a/src/config/option_collation.hpp b/src/config/option_collation.hpp index 1ef016864..20c07419c 100644 --- a/src/config/option_collation.hpp +++ b/src/config/option_collation.hpp @@ -88,6 +88,7 @@ bool annotate_filter_output(const OptionMap& options); boost::optional bamout_request(const OptionMap& options); boost::optional split_bamout_request(const OptionMap& options); +bool full_bamouts_requested(const OptionMap& options); unsigned estimate_max_open_files(const OptionMap& options); diff --git a/src/config/option_parser.cpp b/src/config/option_parser.cpp index 169c60eb0..996107274 100644 --- a/src/config/option_parser.cpp +++ b/src/config/option_parser.cpp @@ -184,13 +184,17 @@ OptionMap parse_options(const int argc, const char** argv) po::value(), "Output realigned BAM files") - ("split-bamout", - po::value(), - "Output split realigned BAM files") + ("split-bamout", + po::value(), + "Output split realigned BAM files") - ("data-profile", - po::value(), - "Output a profile of polymorphisms and errors found in the data") + ("full-bamout", + po::bool_switch()->default_value(false), + "Output all reads when producing realigned bam outputs rather than just variant read minibams") + + ("data-profile", + po::value(), + "Output a profile of polymorphisms and errors found in the data") ; po::options_description transforms("Read transformations"); diff --git a/src/core/calling_components.cpp b/src/core/calling_components.cpp index 19819ec04..099ddb584 100644 --- a/src/core/calling_components.cpp +++ b/src/core/calling_components.cpp @@ -168,6 +168,11 @@ boost::optional GenomeCallingComponents::split_ba return components_.split_bamout; } +bool GenomeCallingComponents::write_full_bamouts() const noexcept +{ + return components_.write_full_bamouts; +} + boost::optional GenomeCallingComponents::data_profile() const { return components_.data_profile; @@ -501,6 +506,7 @@ GenomeCallingComponents::Components::Components(ReferenceGenome&& reference, Rea , filter_request {} , bamout {options::bamout_request(options)} , split_bamout {options::split_bamout_request(options)} +, write_full_bamouts {options::full_bamouts_requested(options)} , data_profile {options::data_profile_request(options)} { drop_unused_samples(this->samples, this->read_manager); diff --git a/src/core/calling_components.hpp b/src/core/calling_components.hpp index f706ef49a..4e1ab8df4 100644 --- a/src/core/calling_components.hpp +++ b/src/core/calling_components.hpp @@ -72,6 +72,7 @@ class GenomeCallingComponents boost::optional filter_request() const; boost::optional bamout() const; boost::optional split_bamout() const; + bool write_full_bamouts() const noexcept; boost::optional data_profile() const; private: @@ -108,7 +109,9 @@ class GenomeCallingComponents bool sites_only; boost::optional legacy; boost::optional filter_request; - boost::optional bamout, split_bamout, data_profile; + boost::optional bamout, split_bamout; + bool write_full_bamouts; + boost::optional data_profile; // Components that require temporary directory during construction appear last to make // exception handling easier. boost::optional temp_directory; diff --git a/src/core/octopus.cpp b/src/core/octopus.cpp index a9e57c478..fa3d704e8 100644 --- a/src/core/octopus.cpp +++ b/src/core/octopus.cpp @@ -1551,9 +1551,11 @@ void run_bam_realign(GenomeCallingComponents& components) if (is_bam_realignment_requested(components)) { if (check_bam_realign(components)) { components.read_manager().close(); + BAMRealigner::Config config {}; + config.copy_hom_ref_reads = components.write_full_bamouts(); if (components.read_manager().paths().size() == 1) { realign(components.read_manager().paths().front(), get_bam_realignment_vcf(components), - *components.bamout(), components.reference()); + *components.bamout(), components.reference(), config); } else { namespace fs = boost::filesystem; const auto bamout_directory = *components.bamout(); @@ -1574,7 +1576,7 @@ void run_bam_realign(GenomeCallingComponents& components) auto bamout_path = bamout_directory; bamout_path /= bamin_path.filename(); if (bamin_path != bamout_path) { - realign(bamin_path, get_bam_realignment_vcf(components), bamout_path, components.reference()); + realign(bamin_path, get_bam_realignment_vcf(components), bamout_path, components.reference(), config); } else { logging::WarningLogger warn_log {}; stream(warn_log) << "Cannot make evidence bam " << bamout_path << " as it is an input bam"; @@ -1586,10 +1588,12 @@ void run_bam_realign(GenomeCallingComponents& components) if (is_split_bam_realignment_requested(components)) { if (check_bam_realign(components)) { components.read_manager().close(); + BAMRealigner::Config config {}; + config.copy_hom_ref_reads = components.write_full_bamouts(); if (components.read_manager().paths().size() == 1) { auto out_paths = get_haplotype_bam_paths(*components.split_bamout(), get_max_ploidy(components)); realign(components.read_manager().paths().front(), get_bam_realignment_vcf(components), - std::move(out_paths), components.reference()); + std::move(out_paths), components.reference(), config); } else { namespace fs = boost::filesystem; const auto bamout_directory = *components.split_bamout(); @@ -1612,7 +1616,7 @@ void run_bam_realign(GenomeCallingComponents& components) bamout_prefix /= bamin_path.filename().stem(); const auto max_ploidy = get_max_called_ploidy(realignment_vcf, bamin_path); auto bamout_paths = get_haplotype_bam_paths(bamout_prefix, max_ploidy); - realign(bamin_path, realignment_vcf, std::move(bamout_paths), components.reference()); + realign(bamin_path, realignment_vcf, std::move(bamout_paths), components.reference(), config); } } }