sdsl/include/sdsl/algorithms_for_string_matching.hpp

Go to the documentation of this file.

/* sdsl - succinct data structures library
    Copyright (C) 2010 Simon Gog

    This program is free software: you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    This program is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with this program.  If not, see http://www.gnu.org/licenses/ .
*/
#ifndef INCLUDED_SDSL_ALGORITHMS_FOR_STRING_MATCHING
#define INCLUDED_SDSL_ALGORITHMS_FOR_STRING_MATCHING

#include "int_vector.hpp"

#include <stdexcept> // for exceptions
#include <iostream>
#include <cassert>
#include <stack>
#include <utility>


namespace sdsl
{

namespace algorithm
{


template<class Csa>
static typename Csa::csa_size_type backward_search(
    const Csa& csa,
    typename Csa::size_type l,
    typename Csa::size_type r,
    typename Csa::char_type c,
    typename Csa::size_type& l_res,
    typename Csa::size_type& r_res
)
{
    assert(l <= r); assert(r < csa.size());
    typename Csa::size_type c_begin = csa.C[csa.char2comp[c]];
    l_res = c_begin + csa.rank_bwt(l, c); // count c in bwt[0..l-1]
    r_res = c_begin + csa.rank_bwt(r+1, c) - 1; // count c in bwt[0..r]
    assert(r_res+1-l_res >= 0);
    return r_res+1-l_res;
}


template<class Csa>
static typename Csa::csa_size_type backward_search(const Csa& csa, typename Csa::size_type l,
        typename Csa::size_type r, typename Csa::pattern_type pat,
        typename Csa::size_type len, typename Csa::size_type& l_res, typename Csa::size_type& r_res)
{
    typename Csa::size_type i = 0;
    l_res = l; r_res = r;
    while (i < len and backward_search(csa, l_res, r_res, pat[len-i-1], l_res, r_res)) {
        ++i;
    }
    return r_res+1-l_res;
}

// TODO: forward search. Include original text???


template<class Csa>
static typename Csa::csa_size_type count(const Csa& csa, typename Csa::pattern_type pat, typename Csa::size_type len)
{
    if (len > csa.size())
        return 0;
    typename Csa::size_type t1,t2; t1=t2=0; // dummy varaiable for the backward_search call
    return backward_search(csa, 0, csa.size()-1, pat, len, t1, t2);
}


template<class Csa, class RandomAccessContainer>
static typename Csa::csa_size_type locate(const Csa&  csa, typename Csa::pattern_type pattern,
        typename Csa::size_type pattern_len, RandomAccessContainer& occ)
{
    typename Csa::size_type occ_begin, occ_end, occs;
    occs = backward_search(csa, 0, csa.size()-1, pattern, pattern_len, occ_begin, occ_end);
    occ.resize(occs);
    for (typename Csa::size_type i=0; i < occs; ++i) {
        occ[i] = csa[occ_begin+i];
    }
    return occs;
}


// Is it cheaper to call T[i] = BWT[iSA[i+1]]??? Additional ranks but H_0 average access
// TODO: extract backward!!! is faster in most cases!
template<class Csa>
static void extract(const Csa& csa, typename Csa::size_type begin, typename Csa::size_type end, unsigned char* text)
{
    assert(end <= csa.size());
    assert(begin <= end);
    for (typename Csa::size_type i=begin, order = csa(begin); i<=end; ++i, order =  csa.psi[order]) {
        uint16_t c_begin = 1, c_end = 257, mid;
        while (c_begin < c_end) {
            mid = (c_begin+c_end)>>1;
            if (csa.C[mid] <= order) {
                c_begin = mid+1;
            } else {
                c_end = mid;
            }
        }
        text[i-begin] = csa.comp2char[c_begin-1];
    }
    if (text[end-begin]!=0)
        text[end-begin+1] = 0; // set terminal character
}


// TODO: use extract with four parameters to implement this method
template<class Csa>
static std::string extract(const Csa& csa, typename Csa::size_type begin, typename Csa::size_type end)
{
    assert(end <= csa.size());
    assert(begin <= end);
    std::string result(end-begin+1,' ');
    for (typename Csa::size_type i=begin, order = csa(begin); i<=end; ++i, order =  csa.psi[order]) {
        uint16_t c_begin = 1, c_end = 257, mid;
        while (c_begin < c_end) {
            mid = (c_begin+c_end)>>1;
            if (csa.C[mid] <= order) {
                c_begin = mid+1;
            } else {
                c_end = mid;
            }
        }
        result[i-begin] = csa.comp2char[c_begin-1];
    }
    return result;
}


template<class Cst>
typename Cst::cst_size_type forward_search(const Cst& cst, typename Cst::node_type& v,
        const typename Cst::size_type d, const typename Cst::char_type c,
        typename Cst::size_type& char_pos)
{
    unsigned char cc = cst.csa.char2comp[c]; // check if c occures in the text of the csa
    if (cc==0 and cc!=c)                                                 //   "    " "    "     "  "    "   "  "   "
        return 0;
    typename Cst::size_type depth_node = cst.depth(v);
    if (d < depth_node) {  //
        char_pos = cst.csa.psi[char_pos];
        if (char_pos < cst.csa.C[cc] or char_pos >= cst.csa.C[cc+1])
            return 0;
        return cst.leaves_in_the_subtree(v);
    } else if (d == depth_node) { //
        v = cst.child(v, c, char_pos);
        if (v == cst.root())
            return 0;
        else
            return cst.leaves_in_the_subtree(v);
    } else {
        throw std::invalid_argument("depth d should be smaller or equal than node depth of v!");
        return 0;
    }
}


template<class Cst>
typename Cst::cst_size_type forward_search(const Cst& cst, typename Cst::node_type& v,
        typename Cst::size_type d, typename Cst::pattern_type pat, typename Cst::size_type len,
        typename Cst::size_type& char_pos
                                          )
{
    if (len==0)
        return cst.leaves_in_the_subtree(v);
    typename Cst::size_type i=0, size=0;
    while (i < len and (size=forward_search(cst, v, d, pat[i], char_pos))) {
        ++d;
        ++i;
    }
    return size;
}


template<class Cst>
typename Cst::cst_size_type count(const Cst& cst, typename Cst::pattern_type pat, typename Cst::size_type len)
{
    typename Cst::node_type v = cst.root();
    typename Cst::size_type char_pos = 0;
    return forward_search(cst, v, 0, pat, len, char_pos);
}


template<class Cst, class RandomAccessContainer>
typename Cst::cst_size_type locate(const Cst& cst, typename Cst::pattern_type pat,
                                   typename Cst::size_type len, RandomAccessContainer& occ)
{
    typedef typename Cst::size_type size_type;
    typename Cst::node_type v = cst.root();
    size_type char_pos = 0;
    typename Cst::cst_size_type occs = forward_search(cst, v, 0, pat, len, char_pos);
    occ.resize(occs);
    if (occs == 0)
        return 0; // because v equals cst.root()
    size_type left      = cst.lb(v);   // get range in the
    size_type right = cst.rb(v);
    for (size_type i=left; i <= right; ++i)
        occ[i-left] = cst.csa[i];
    return occs;
}


template<class Cst>
void extract(const Cst& cst, const typename Cst::node_type& v, unsigned char* text)
{
    if (v == cst.root()) {
        text[0] = 0;
        return;
    }
    // first get the suffix array entry of the leftmost leaf in the subtree rooted at v
    typename Cst::size_type begin = cst.csa[cst.lb(v)];
    // then call the extract method on the compressed suffix array
    extract(cst.csa, begin, begin + cst.depth(v) - 1, text);
}


template<class Cst>
std::string extract(const Cst& cst, const typename Cst::node_type& v)
{
    if (v==cst.root()) {
        return std::string();
    }
    // first get the suffix array entry of the leftmost leaf in the subtree rooted at v
    typename Cst::size_type begin = cst.csa[cst.lb(v)];
    // then call the extract method on the compressed suffix array
    return extract(cst.csa, begin, begin + cst.depth(v) - 1);
}


/*
template<class Cst>
typename Cst::size_type count(
                        const Cst &cst,
                        typename Cst::pattern_type pattern,
                    typename Cst::size_type pattern_len){
        if(pattern_len==0){
                return 0;
        }
        typedef typename Cst::size_type size_type;
        typedef typename Cst::node_type node_type;
        node_type node = cst.root();
        for(size_type i=0, char_pos=0; cst.depth(node) < pattern_len; ++i){
                node_type newnode = cst.child(node, (unsigned char)pattern[cst.depth(node)], char_pos);
                if( newnode == cst.root() )// root node, no match found
                        return 0;
                // else the first character of the newnode matches the pattern at position depth(node)
                for(size_type j=cst.depth(node)+1; j < cst.depth(newnode) and j < pattern_len; ++j){
                        char_pos = cst.csa.psi[char_pos];
                        size_type cc = cst.csa.char2comp[pattern[j]];
                        if(char_pos < cst.csa.C[cc] or char_pos >= cst.csa.C[cc+1] )
                                return 0;
                }
                node = newnode;
        }
        return cst.leaves_in_the_subtree(node);
}
*/
/*
template<class Cst, class RandomAccessContainer>
typename Cst::size_type locate(
                const Cst &cst,
                typename Cst::pattern_type pattern,
                typename Cst::size_type pattern_len,
                RandomAccessContainer &occ){
        occ.resize(0);
        typedef typename Cst::size_type size_type;
        typedef typename Cst::node_type node_type;
        node_type node = cst.root();
        for(size_type i=0, char_pos=0; cst.depth(node) < pattern_len; ++i){
                node_type newnode = cst.child(node, (unsigned char)pattern[cst.depth(node)], char_pos);
                if( newnode == cst.root() )// root node, no match found
                        return 0;
                // else the first character of the newnode matches the pattern at position depth(node)
                for(size_type j=cst.depth(node)+1; j < cst.depth(newnode) and j < pattern_len; ++j){
                        char_pos = cst.csa.psi[char_pos];
                        size_type cc = cst.csa.char2comp[pattern[j]];
                        if(char_pos < cst.csa.C[cc] or char_pos >= cst.csa.C[cc+1] )
                                return 0;
                }
                node = newnode;
        }
        size_type occs = cst.leaves_in_the_subtree(node);
        occ.resize(occs);
        size_type left  = cst.leftmost_suffix_array_index_in_the_subtree(node);
        size_type right = cst.rightmost_suffix_array_index_in_the_subtree(node);
        for(size_type i=left; i <= right; ++i)
                occ[i-left] = cst.csa[i];
        return occs;
}
*/


/*
template<class Csa, class Lcp, class Bp_support>
typename cst_sct<Csa, Lcp, Bp_support>::size_type count(
                        const cst_sct<Csa, Lcp, Bp_support> &cst,
                        typename cst_sct<Csa, Lcp, Bp_support>::pattern_type pattern,
                    typename cst_sct<Csa, Lcp, Bp_support>::size_type pattern_len){
        if(pattern_len==0){
                return 0;
        }
        typedef typename cst_sct<Csa, Lcp, Bp_support>::size_type size_type;
        typedef typename cst_sct<Csa, Lcp, Bp_support>::node_type node_type;
        node_type node = cst.root();
        for(size_type i=0, char_pos=0; node.l < pattern_len; ++i){
                node_type newnode = cst.child(node, (unsigned char)pattern[node.l], char_pos);
                if( newnode.l == 0 )// root node, no match found
                        return 0;
                // else the first character of the newnode matches the pattern at position node.l
                for(size_type j=node.l+1; j < newnode.l and j< pattern_len; ++j){
                        char_pos = cst.csa.psi[char_pos];
                        size_type cc = cst.csa.char2comp[pattern[j]];
                        if(char_pos < cst.csa.C[cc] or char_pos >= cst.csa.C[cc+1] )
                                return 0;
                }
                node = newnode;
        }
        return cst.leaves_in_the_subtree(node);
}
*/

} // end namespace algorithm

} // end namespace sdsl

#endif