structure_file/input.hpp

Go to the documentation of this file.

// SPDX-FileCopyrightText: 2006-2025 Knut Reinert & Freie Universität Berlin
// SPDX-FileCopyrightText: 2016-2025 Knut Reinert & MPI für molekulare Genetik
// SPDX-License-Identifier: BSD-3-Clause

 
#pragma once
 
#include <cassert>
#include <filesystem>
#include <fstream>
#include <limits>
#include <optional>
#include <string>
#include <type_traits>
#include <utility>
#include <variant>
#include <vector>
 
#include <seqan3/alphabet/adaptation/char.hpp>
#include <seqan3/alphabet/aminoacid/aa27.hpp>
#include <seqan3/alphabet/nucleotide/rna15.hpp>
#include <seqan3/alphabet/nucleotide/rna5.hpp>
#include <seqan3/alphabet/structure/dssp9.hpp>
#include <seqan3/alphabet/structure/structured_aa.hpp>
#include <seqan3/io/detail/in_file_iterator.hpp>
#include <seqan3/io/detail/misc_input.hpp>
#include <seqan3/io/detail/record.hpp>
#include <seqan3/io/exception.hpp>
#include <seqan3/io/stream/concept.hpp>
#include <seqan3/io/structure_file/format_vienna.hpp>
#include <seqan3/io/structure_file/input_format_concept.hpp>
#include <seqan3/io/structure_file/input_options.hpp>
#include <seqan3/io/structure_file/record.hpp>
#include <seqan3/utility/container/concept.hpp>
#include <seqan3/utility/type_list/traits.hpp>
 
namespace seqan3
{
// ----------------------------------------------------------------------------
// structure_file_input_traits
// ----------------------------------------------------------------------------

template <typename t>
concept structure_file_input_traits = requires (t v) {
    // TODO(joergi-w) The expensive concept checks are currently omitted. Check again when compiler has improved.
    // sequence
    requires writable_alphabet<typename t::seq_alphabet>;
    requires writable_alphabet<typename t::seq_legal_alphabet>;
    requires explicitly_convertible_to<typename t::seq_legal_alphabet, typename t::seq_alphabet>;
    requires sequence_container<typename t::template seq_container<typename t::seq_alphabet>>;
 
    // id
    requires writable_alphabet<typename t::id_alphabet>;
    requires sequence_container<typename t::template id_container<typename t::id_alphabet>>;
 
    // bpp
    requires std::is_floating_point_v<typename t::bpp_prob>;
    requires std::numeric_limits<typename t::bpp_partner>::is_integer;
 
    //    requires container // TODO check Associative container Concept when implemented
    //        <typename t::template bpp_queue
    //            <typename t::template bpp_item
    //                <typename t::bpp_prob, typename t::bpp_partner>>>
    //        && requires(typename t::template bpp_queue // TODO maybe implement also a version that allows emplace_back
    //            <typename t::template bpp_item
    //                <typename t::bpp_prob, typename t::bpp_partner>> value) { value.emplace(1.0, 1); };
    //    requires sequence_container
    //        <typename t::template bpp_container
    //            <typename t::template bpp_queue
    //                 <typename t::template bpp_item
    //                      <typename t::bpp_prob, typename t::bpp_partner>>>>;
 
    // structure
    requires std::is_same_v<typename t::structure_alphabet, dssp9> // TODO(joergi-w) add aa_structure_concept
                 || rna_structure_alphabet<typename t::structure_alphabet>;
    requires sequence_container<typename t::template structure_container<typename t::structure_alphabet>>;
 
    // structured sequence: tuple composites of seq and structure
    requires std::is_base_of_v<
        alphabet_tuple_base<
            typename t::template structured_seq_alphabet<typename t::seq_alphabet, typename t::structure_alphabet>,
            typename t::seq_alphabet,
            typename t::structure_alphabet>,
        typename t::template structured_seq_alphabet<typename t::seq_alphabet, typename t::structure_alphabet>>;
    //    requires sequence_container
    //        <typename t::template structured_seq_container
    //            <typename t::template structured_seq_alphabet
    //                <typename t::seq_alphabet, typename t::structure_alphabet>>>;
 
    // energy: std::optional of floating point number
    requires std::is_floating_point_v<typename t::energy_type::value_type>;
 
    // reactivity [error]
    requires std::is_floating_point_v<typename t::react_type>;
    requires sequence_container<typename t::template react_container<typename t::react_type>>;
 
    // comment
    requires writable_alphabet<typename t::comment_alphabet>;
    requires sequence_container<typename t::template comment_container<typename t::comment_alphabet>>;
 
    // offset
    requires std::numeric_limits<typename t::offset_type>::is_integer;
};
 
// ----------------------------------------------------------------------------
// structure_file_input_default_traits
// ----------------------------------------------------------------------------

struct structure_file_input_default_traits_rna
{
 
    // sequence

    using seq_alphabet = rna5;

    using seq_legal_alphabet = rna15;

    template <typename _seq_alphabet>
    using seq_container = std::vector<_seq_alphabet>;
 
    // id

    using id_alphabet = char;

    template <typename _id_alphabet>
    using id_container = std::basic_string<_id_alphabet>;
 
    // base pair probability structure

    using bpp_prob = double;

    using bpp_partner = size_t;

    template <typename _bpp_prob, typename _bpp_partner>
    using bpp_item = std::pair<_bpp_prob, _bpp_partner>;

    template <typename _bpp_item>
    using bpp_queue = std::set<_bpp_item>;

    template <typename _bpp_queue>
    using bpp_container = std::vector<_bpp_queue>;
 
    // fixed structure

    using structure_alphabet = wuss51;

    template <typename _structure_alphabet>
    using structure_container = std::vector<_structure_alphabet>;
 
    // combined sequence and structure

    template <typename _seq_alphabet, typename _structure_alphabet>
    using structured_seq_alphabet = structured_rna<_seq_alphabet, _structure_alphabet>;

    template <typename _structured_seq_alphabet>
    using structured_seq_container = std::vector<_structured_seq_alphabet>;
 
    // energy

    using energy_type = std::optional<double>;
 
    // reactivity [error]

    using react_type = double;

    template <typename _react_type>
    using react_container = std::vector<_react_type>;
 
    // comment

    using comment_alphabet = char;

    template <typename _comment_alphabet>
    using comment_container = std::basic_string<_comment_alphabet>;
 
    // offset

    using offset_type = size_t;
};

struct structure_file_input_default_traits_aa : structure_file_input_default_traits_rna
{

    using seq_alphabet = aa27;
    using seq_legal_alphabet = aa27;
    using structure_alphabet = dssp9;
    template <typename _seq_alphabet, typename _structure_alphabet>
    using structured_seq_alphabet = structured_aa<_seq_alphabet, _structure_alphabet>;
};
 
// ----------------------------------------------------------------------------
// structure_file_input
// ----------------------------------------------------------------------------

template <structure_file_input_traits traits_type_ = structure_file_input_default_traits_rna,
          detail::fields_specialisation selected_field_ids_ = fields<field::seq, field::id, field::structure>,
          detail::type_list_of_structure_file_input_formats valid_formats_ = type_list<format_vienna>>
class structure_file_input
{
public:
    using traits_type = traits_type_;
    using selected_field_ids = selected_field_ids_;
    using valid_formats = valid_formats_;
    using stream_char_type = char;

    using field_ids = fields<field::seq,
                             field::id,
                             field::bpp,
                             field::structure,
                             field::structured_seq,
                             field::energy,
                             field::react,
                             field::react_err,
                             field::comment,
                             field::offset>;
 
    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 structure files, please refer to the documentation "
        "of structure_file_input::field_ids for the accepted values.");
 
    static_assert(
        []() constexpr
        {
            return !(selected_field_ids::contains(field::structured_seq)
                     && (selected_field_ids::contains(field::seq) || (selected_field_ids::contains(field::structure))));
        }(),
        "You may not select field::structured_seq and either of field::seq and field::structure "
        "at the same time.");

    using seq_type = typename traits_type::template seq_container<typename traits_type::seq_alphabet>;
    using id_type = typename traits_type::template id_container<typename traits_type::id_alphabet>;
    using bpp_type = typename traits_type::template bpp_container<typename traits_type::template bpp_queue<
        typename traits_type::template bpp_item<typename traits_type::bpp_prob, typename traits_type::bpp_partner>>>;
    using structure_type = typename traits_type::template structure_container<typename traits_type::structure_alphabet>;
    using structured_seq_type = typename traits_type::template structured_seq_container<
        typename traits_type::template structured_seq_alphabet<typename traits_type::seq_alphabet,
                                                               typename traits_type::structure_alphabet>>;
    using energy_type = typename traits_type::energy_type;
    using react_type = typename traits_type::template react_container<typename traits_type::react_type>;
    using comment_type = typename traits_type::template comment_container<typename traits_type::comment_alphabet>;
    using offset_type = typename traits_type::offset_type;

    using field_types = type_list<seq_type,
                                  id_type,
                                  bpp_type,
                                  structure_type,
                                  structured_seq_type,
                                  energy_type,
                                  react_type,
                                  react_type,
                                  comment_type,
                                  offset_type>;

    using record_type = structure_record<detail::select_types_with_ids_t<field_types, field_ids, selected_field_ids>,
                                         selected_field_ids>;

    using value_type = record_type;
    using reference = record_type &;
    using const_reference = void;
    using size_type = size_t;
    using difference_type = std::make_signed_t<size_t>;
    using iterator = detail::in_file_iterator<structure_file_input>;
    using const_iterator = void;
    using sentinel = std::default_sentinel_t;

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

    structure_file_input(std::filesystem::path filename,
                         selected_field_ids const & SEQAN3_DOXYGEN_ONLY(fields_tag) = selected_field_ids{}) :
        primary_stream{new std::ifstream{}, stream_deleter_default}
    {
        primary_stream->rdbuf()->pubsetbuf(stream_buffer.data(), stream_buffer.size());
        static_cast<std::basic_ifstream<char> *>(primary_stream.get())
            ->open(filename, std::ios_base::in | std::ios::binary);
 
        if (!primary_stream->good())
            throw file_open_error{"Could not open file " + filename.string() + " for reading."};
 
        // possibly add intermediate decompression stream
        secondary_stream = detail::make_secondary_istream(*primary_stream, filename);
 
        // initialise format handler
        detail::set_format(format, filename);
    }

    template <input_stream stream_t, structure_file_input_format file_format>
        requires std::same_as<typename std::remove_reference_t<stream_t>::char_type, char>
    structure_file_input(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},
        format{detail::structure_file_input_format_exposer<file_format>{}}
    {
        static_assert(list_traits::contains<file_format, valid_formats>,
                      "You selected a format that is not in the valid_formats of this file.");
 
        // possibly add intermediate decompression stream
        secondary_stream = detail::make_secondary_istream(*primary_stream);
    }

    template <input_stream stream_t, structure_file_input_format file_format>
        requires std::same_as<typename std::remove_reference_t<stream_t>::char_type, char>
    structure_file_input(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},
        format{detail::structure_file_input_format_exposer<file_format>{}}
    {
        static_assert(list_traits::contains<file_format, valid_formats>,
                      "You selected a format that is not in the valid_formats of this file.");
 
        // possibly add intermediate compression stream
        secondary_stream = detail::make_secondary_istream(*primary_stream);
    }


    iterator begin()
    {
        // buffer first record
        if (!first_record_was_read)
        {
            read_next_record();
            first_record_was_read = true;
        }
 
        return {*this};
    }

    sentinel end() noexcept
    {
        return {};
    }

    reference front() noexcept
    {
        return *begin();
    }


    structure_file_input_options<typename traits_type::seq_legal_alphabet,
                                 selected_field_ids::contains(field::structured_seq)>
        options;
 
protected:

    record_type record_buffer;
    std::vector<char> stream_buffer{std::vector<char>(1'000'000)};
    std::streampos position_buffer{};

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

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

    bool first_record_was_read{false};
    bool at_end{false};

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

    void read_next_record()
    {
        // clear the record
        record_buffer.clear();
 
        // store the current position in the position buffer
        position_buffer = secondary_stream->tellg();
 
        // at end if we could not read further
        if ((std::istreambuf_iterator<stream_char_type>{*secondary_stream}
             == std::istreambuf_iterator<stream_char_type>{}))
        {
            at_end = true;
            return;
        }
 
        assert(!format.valueless_by_exception());
        std::visit(
            [&](auto & f)
            {
                // read new record
                if constexpr (selected_field_ids::contains(field::structured_seq))
                {
                    static_assert(!selected_field_ids::contains(field::structure),
                                  "You may not select field::structured_seq and field::structure at the same time.");
                    static_assert(!selected_field_ids::contains(field::seq),
                                  "You may not select field::structured_seq and field::seq at the same time.");
                    f.read_structure_record(*secondary_stream,
                                            options,
                                            detail::get_or_ignore<field::structured_seq>(record_buffer), // seq
                                            detail::get_or_ignore<field::id>(record_buffer),
                                            detail::get_or_ignore<field::bpp>(record_buffer),
                                            detail::get_or_ignore<field::structured_seq>(record_buffer), // structure
                                            detail::get_or_ignore<field::energy>(record_buffer),
                                            detail::get_or_ignore<field::react>(record_buffer),
                                            detail::get_or_ignore<field::react_err>(record_buffer),
                                            detail::get_or_ignore<field::comment>(record_buffer),
                                            detail::get_or_ignore<field::offset>(record_buffer));
                }
                else
                {
                    f.read_structure_record(*secondary_stream,
                                            options,
                                            detail::get_or_ignore<field::seq>(record_buffer),
                                            detail::get_or_ignore<field::id>(record_buffer),
                                            detail::get_or_ignore<field::bpp>(record_buffer),
                                            detail::get_or_ignore<field::structure>(record_buffer),
                                            detail::get_or_ignore<field::energy>(record_buffer),
                                            detail::get_or_ignore<field::react>(record_buffer),
                                            detail::get_or_ignore<field::react_err>(record_buffer),
                                            detail::get_or_ignore<field::comment>(record_buffer),
                                            detail::get_or_ignore<field::offset>(record_buffer));
                }
            },
            format);
    }

    friend iterator;
};


template <input_stream stream_type,
          structure_file_input_format file_format,
          detail::fields_specialisation selected_field_ids>
structure_file_input(stream_type && stream, file_format const &, selected_field_ids const &)
    -> structure_file_input<typename structure_file_input<>::traits_type, // actually use the default
                            selected_field_ids,
                            type_list<file_format>>;

template <input_stream stream_type,
          structure_file_input_format file_format,
          detail::fields_specialisation selected_field_ids>
structure_file_input(stream_type & stream, file_format const &, selected_field_ids const &)
    -> structure_file_input<typename structure_file_input<>::traits_type, // actually use the default
                            selected_field_ids,
                            type_list<file_format>>;
 
} // namespace seqan3