detail.hpp

Go to the documentation of this file.

// -----------------------------------------------------------------------------------------------------
// Copyright (c) 2006-2022, Knut Reinert & Freie Universität Berlin
// Copyright (c) 2016-2022, 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 <map>
#include <ranges>
#include <stack>
 
#include <seqan3/alphabet/structure/concept.hpp>
#include <seqan3/io/exception.hpp>
 
namespace seqan3::detail
{
template <typename structure_alph_type, typename bpp_type, std::ranges::range structure_type>
inline void bpp_from_rna_structure(bpp_type & bpp, structure_type const & structure, double weight = 1.)
{
    if constexpr (!rna_structure_alphabet<structure_alph_type>)
        throw parse_error{"Cannot create base pair probabilities from a structure that is not RNA structure."};
 
    bpp.clear();
    if constexpr (std::ranges::sized_range<structure_type>)
        bpp.reserve(size(structure));
 
    std::stack<size_t> brackets[max_pseudoknot_depth<structure_alph_type>];
    size_t pos = 0ul;
    for (structure_alph_type symbol : structure)
    {
        bpp.push_back({});
        uint8_t const id = pseudoknot_id(symbol).value_or(0);
 
        if (symbol.is_pair_open())
        {
            brackets[id].push(pos);
        }
        else if (symbol.is_pair_close())
        {
            if (!brackets[id].empty())
            {
                bpp[pos].emplace(weight, brackets[id].top());
                bpp[brackets[id].top()].emplace(weight, pos);
                brackets[id].pop();
            }
            else
            {
                throw parse_error{std::string{"Invalid bracket notation: Unpaired closing bracket at position "}
                                  + std::to_string(pos) + "."};
            };
        }
        // no actions for unpaired
        ++pos;
    }
    for (uint8_t id = 0u; id < max_pseudoknot_depth<structure_alph_type>; ++id)
    {
        if (!brackets[id].empty())
        {
            throw parse_error{std::string{"Invalid bracket notation: Unpaired opening bracket at position "}
                              + std::to_string(brackets[id].top()) + "."};
        }
    }
}
 
} // namespace seqan3::detail