output.hpp

Go to the documentation of this file.

// -----------------------------------------------------------------------------------------------------
// Copyright (c) 2006-2019, Knut Reinert & Freie Universität Berlin
// Copyright (c) 2016-2019, Knut Reinert & MPI für molekulare Genetik
// This file may be used, modified and/or redistributed under the terms of the 3-clause BSD-License
// shipped with this file and also available at: https://github.com/seqan/seqan3/blob/master/LICENSE.md
// -----------------------------------------------------------------------------------------------------

#pragma once

#include <cassert>
#include <fstream>
#include <string>
#include <variant>
#include <vector>

#include <range/v3/algorithm/equal.hpp>
#include <range/v3/view/zip.hpp>

#include <seqan3/core/type_traits/basic.hpp>
#include <seqan3/core/type_traits/template_inspection.hpp>
#include <seqan3/core/concept/tuple.hpp>
#include <seqan3/io/stream/concept.hpp>
#include <seqan3/io/exception.hpp>
#include <seqan3/std/filesystem>
#include <seqan3/io/record.hpp>
#include <seqan3/io/detail/misc.hpp>
#include <seqan3/io/detail/misc_output.hpp>
#include <seqan3/io/detail/out_file_iterator.hpp>
#include <seqan3/io/detail/record.hpp>
#include <seqan3/io/sequence_file/format_embl.hpp>
#include <seqan3/io/sequence_file/format_fasta.hpp>
#include <seqan3/io/sequence_file/format_fastq.hpp>
#include <seqan3/io/sequence_file/format_genbank.hpp>
#include <seqan3/io/sequence_file/format_sam.hpp>
#include <seqan3/io/sequence_file/output_format_concept.hpp>
#include <seqan3/io/sequence_file/output_options.hpp>
#include <seqan3/range/view/convert.hpp>
#include <seqan3/range/view/get.hpp>
#include <seqan3/std/ranges>

namespace seqan3
{

// ----------------------------------------------------------------------------
// sequence_file_output
// ----------------------------------------------------------------------------

template <detail::Fields selected_field_ids_ = fields<field::SEQ, field::ID, field::QUAL>,
          detail::TypeListOfSequenceFileOutputFormats valid_formats_ =
              type_list<format_embl, format_fasta, format_fastq, format_genbank, format_sam>,
          Char stream_char_type_ = char>
class sequence_file_output
{
public:
    using selected_field_ids    = selected_field_ids_;
    using valid_formats         = valid_formats_;
    using stream_char_type      = stream_char_type_;

    using field_ids            = fields<field::SEQ, field::ID, field::QUAL, field::SEQ_QUAL>;

    static_assert([] () constexpr
                  {
                      for (field f : selected_field_ids::as_array)
                          if (!field_ids::contains(f))
                              return false;
                      return true;
                  }(),
                  "You selected a field that is not valid for sequence files, please refer to the documentation "
                  "of sequence_file_output::field_ids for the accepted values.");

    static_assert([] () constexpr
                  {
                      return !(selected_field_ids::contains(field::SEQ_QUAL) &&
                               (selected_field_ids::contains(field::SEQ) ||
                               (selected_field_ids::contains(field::QUAL))));
                  }(),
                  "You may not select field::SEQ_QUAL and either of field::SEQ and field::QUAL at the same time.");

    using value_type        = void;
    using reference         = void;
    using const_reference   = void;
    using size_type         = void;
    using difference_type   = std::ptrdiff_t;
    using iterator          = detail::out_file_iterator<sequence_file_output>;
    using const_iterator    = void;
    using sentinel          = std::ranges::default_sentinel_t;

    sequence_file_output() = delete;
    sequence_file_output(sequence_file_output const &) = delete;
    sequence_file_output & operator=(sequence_file_output const &) = delete;
    sequence_file_output(sequence_file_output &&) = default;
    sequence_file_output & operator=(sequence_file_output &&) = default;
    ~sequence_file_output() = default;

    sequence_file_output(std::filesystem::path filename,
                         selected_field_ids const & SEQAN3_DOXYGEN_ONLY(fields_tag) = selected_field_ids{}) :
        primary_stream{new std::ofstream{filename, std::ios_base::out | std::ios::binary}, stream_deleter_default}
    {
        if (!primary_stream->good())
            throw file_open_error{"Could not open file " + filename.string() + " for writing."};

        // possibly add intermediate compression stream
        secondary_stream = detail::make_secondary_ostream(*primary_stream, filename);

        // initialise format handler or throw if format is not found
        detail::set_format(format, filename);
    }

    template <OStream2 stream_t,
              SequenceFileOutputFormat file_format>
    sequence_file_output(stream_t                 & stream,
                         file_format        const & SEQAN3_DOXYGEN_ONLY(format_tag),
                         selected_field_ids const & SEQAN3_DOXYGEN_ONLY(fields_tag) = selected_field_ids{}) :
        primary_stream{&stream, stream_deleter_noop},
        secondary_stream{&stream, stream_deleter_noop},
        format{detail::sequence_file_output_format<file_format>{}}
    {
        static_assert(meta::in<valid_formats, file_format>::value,
                      "You selected a format that is not in the valid_formats of this file.");
    }

    template <OStream2 stream_t,
              SequenceFileOutputFormat file_format>
    sequence_file_output(stream_t                && stream,
                         file_format        const & SEQAN3_DOXYGEN_ONLY(format_tag),
                         selected_field_ids const & SEQAN3_DOXYGEN_ONLY(fields_tag) = selected_field_ids{}) :
        primary_stream{new stream_t{std::move(stream)}, stream_deleter_default},
        secondary_stream{&*primary_stream, stream_deleter_noop},
        format{detail::sequence_file_output_format<file_format>{}}
    {
        static_assert(meta::in<valid_formats, file_format>::value,
                      "You selected a format that is not in the valid_formats of this file.");
    }

    iterator begin() noexcept
    {
        return {*this};
    }

    sentinel end() noexcept
    {
        return {};
    }

    template <typename record_t>
    void push_back(record_t && r)
        requires TupleLike<record_t> &&
                 requires { requires detail::is_type_specialisation_of_v<remove_cvref_t<record_t>, record>; }
    {
        write_record(detail::get_or_ignore<field::SEQ>(r),
                     detail::get_or_ignore<field::ID>(r),
                     detail::get_or_ignore<field::QUAL>(r),
                     detail::get_or_ignore<field::SEQ_QUAL>(r));

    }

    template <typename tuple_t>
    void push_back(tuple_t && t)
        requires TupleLike<tuple_t>
    {
        // index_of might return npos, but this will be handled well by get_or_ignore (and just return ignore)
        write_record(detail::get_or_ignore<selected_field_ids::index_of(field::SEQ)>(t),
                     detail::get_or_ignore<selected_field_ids::index_of(field::ID)>(t),
                     detail::get_or_ignore<selected_field_ids::index_of(field::QUAL)>(t),
                     detail::get_or_ignore<selected_field_ids::index_of(field::SEQ_QUAL)>(t));
    }

    template <typename arg_t, typename ... arg_types>
    void emplace_back(arg_t && arg, arg_types && ... args)
    {
        push_back(std::tie(arg, args...));
    }

    template <std::ranges::InputRange rng_t>
    sequence_file_output & operator=(rng_t && range)
        requires TupleLike<reference_t<rng_t>>
    {
        for (auto && record : range)
            push_back(std::forward<decltype(record)>(record));
        return *this;
    }

    template <std::ranges::InputRange rng_t>
    friend sequence_file_output & operator|(rng_t && range, sequence_file_output & f)
        requires TupleLike<reference_t<rng_t>>
    {
        f = range;
        return f;
    }

    template <std::ranges::InputRange rng_t>
    friend sequence_file_output operator|(rng_t && range, sequence_file_output && f)
        requires TupleLike<reference_t<rng_t>>
    {
    #if defined(__GNUC__) && (__GNUC__ == 9) // an unreported build problem of GCC9
        for (auto && record : range)
            f.push_back(std::forward<decltype(record)>(record));
    #else // ^^^ workaround | regular solution ↓↓↓
        f = range;
    #endif
        return std::move(f);
    }

    template <typename typelist, typename field_ids>
    sequence_file_output & operator=(record<typelist, field_ids> const & r)
    {
        write_columns(detail::range_wrap_ignore(detail::get_or_ignore<field::SEQ>(r)),
                      detail::range_wrap_ignore(detail::get_or_ignore<field::ID>(r)),
                      detail::range_wrap_ignore(detail::get_or_ignore<field::QUAL>(r)),
                      detail::range_wrap_ignore(detail::get_or_ignore<field::SEQ_QUAL>(r)));
        return *this;
    }

    template <typename ... arg_types>
    sequence_file_output & operator=(std::tuple<arg_types...> const & t)
    {
        // index_of might return npos, but this will be handled well by get_or_ignore (and just return ignore)
        write_columns(
            detail::range_wrap_ignore(detail::get_or_ignore<selected_field_ids::index_of(field::SEQ)>(t)),
            detail::range_wrap_ignore(detail::get_or_ignore<selected_field_ids::index_of(field::ID)>(t)),
            detail::range_wrap_ignore(detail::get_or_ignore<selected_field_ids::index_of(field::QUAL)>(t)),
            detail::range_wrap_ignore(detail::get_or_ignore<selected_field_ids::index_of(field::SEQ_QUAL)>(t)));
        return *this;
    }

    sequence_file_output_options options;

    std::basic_ostream<stream_char_type> & get_stream()
    {
        return *secondary_stream;
    }
protected:

    using stream_ptr_t = std::unique_ptr<std::basic_ostream<stream_char_type>,
                                         std::function<void(std::basic_ostream<stream_char_type>*)>>;
    static void stream_deleter_noop(std::basic_ostream<stream_char_type> *) {}
    static void stream_deleter_default(std::basic_ostream<stream_char_type> * ptr) { delete ptr; }

    stream_ptr_t primary_stream{nullptr, stream_deleter_noop};
    stream_ptr_t secondary_stream{nullptr, stream_deleter_noop};

    using format_type = typename detail::variant_from_tags<valid_formats, detail::sequence_file_output_format>::type;
    format_type format;

    template <typename seq_t, typename id_t, typename qual_t, typename seq_qual_t>
    void write_record(seq_t && seq, id_t && id, qual_t && qual, seq_qual_t && seq_qual)
    {
        static_assert(detail::decays_to_ignore_v<seq_qual_t> ||
                      (detail::decays_to_ignore_v<seq_t> && detail::decays_to_ignore_v<qual_t>),
                      "You may not select field::SEQ_QUAL and either of field::SEQ and field::QUAL at the same time.");

        if constexpr (!detail::decays_to_ignore_v<seq_qual_t>)
            static_assert(detail::is_type_specialisation_of_v<value_type_t<seq_qual_t>, qualified>,
                          "The SEQ_QUAL field must contain a range over the seqan3::qualified alphabet.");

        assert(!format.valueless_by_exception());
        std::visit([&] (auto & f)
        {
            if constexpr (!detail::decays_to_ignore_v<seq_qual_t>)
            {
                f.write(*secondary_stream,
                        options,
                        seq_qual | view::get<0>,
                        id,
                        seq_qual | view::get<1>);
            }
            else
            {
                f.write(*secondary_stream,
                        options,
                        seq,
                        id,
                        qual);
            }
        }, format);
    }

    template <std::ranges::InputRange seqs_t,
              std::ranges::InputRange ids_t,
              std::ranges::InputRange quals_t,
              std::ranges::InputRange seq_quals_t>
    void write_columns(seqs_t       && seqs,
                       ids_t        && ids,
                       quals_t      && quals,
                       seq_quals_t  && seq_quals)
    {
        static_assert(!(detail::decays_to_ignore_v<reference_t<seqs_t>> &&
                        detail::decays_to_ignore_v<reference_t<ids_t>> &&
                        detail::decays_to_ignore_v<reference_t<quals_t>> &&
                        detail::decays_to_ignore_v<reference_t<seq_quals_t>>),
                      "At least one of the columns must not be set to std::ignore.");

        static_assert(detail::decays_to_ignore_v<reference_t<seq_quals_t>> ||
                      (detail::decays_to_ignore_v<reference_t<seqs_t>> &&
                       detail::decays_to_ignore_v<reference_t<quals_t>>),
                      "You may not select field::SEQ_QUAL and either of field::SEQ and field::QUAL at the same time.");

        if constexpr (!detail::decays_to_ignore_v<reference_t<seq_quals_t>>)
            static_assert(detail::is_type_specialisation_of_v<value_type_t<reference_t<seq_quals_t>>, qualified>,
                          "The SEQ_QUAL field must contain a range over the seqan3::qualified alphabet.");

        assert(!format.valueless_by_exception());
        std::visit([&] (auto & f)
        {
            if constexpr (!detail::decays_to_ignore_v<reference_t<seq_quals_t>>)
            {
                auto zipped = std::view::zip(seq_quals, ids);

                for (auto && v : zipped)
                    f.write(*secondary_stream,
                            options,
                            std::get<0>(v) | view::get<0>,
                            std::get<1>(v),
                            std::get<0>(v) | view::get<1>);
            }
            else
            {
                auto zipped = std::view::zip(seqs, ids, quals);

                for (auto && v : zipped)
                    f.write(*secondary_stream, options, std::get<0>(v), std::get<1>(v), std::get<2>(v));
            }
        }, format);
    }

    friend iterator;
};

template <OStream2                 stream_t,
          SequenceFileOutputFormat file_format,
          detail::Fields           selected_field_ids>
sequence_file_output(stream_t &&,
                     file_format const &,
                     selected_field_ids const &)
    -> sequence_file_output<selected_field_ids,
                            type_list<file_format>,
                            typename std::remove_reference_t<stream_t>::char_type>;

template <OStream2                 stream_t,
          SequenceFileOutputFormat file_format,
          detail::Fields           selected_field_ids>
sequence_file_output(stream_t &,
                     file_format const &,
                     selected_field_ids const &)
    -> sequence_file_output<selected_field_ids,
                            type_list<file_format>,
                            typename std::remove_reference_t<stream_t>::char_type>;
} // namespace seqan3