source: trunk/SECEDIT/SEC_helix.cxx

// =============================================================== //
//                                                                 //
//   File      : SEC_helix.cxx                                     //
//   Purpose   : helix related things                              //
//                                                                 //
//   Coded by Ralf Westram (coder@reallysoft.de) in August 2007    //
//   Institute of Microbiology (Technical University Munich)       //
//   http://www.arb-home.de/                                       //
//                                                                 //
// =============================================================== //

#include "SEC_root.hxx"
#include "SEC_helix.hxx"
#include <arbdbt.h>
#include <arb_strbuf.h>

using namespace std;

const size_t *SEC_root::getHelixPositions(const char *helixNr) const {
    sec_assert(helixNr);
    const BI_helix *helix  = get_helixDef();
    long            start1 = helix->first_position(helixNr);

    if (start1 == -1) return NULp;

    size_t        start2 = helix->last_position(helixNr);
    static size_t pos[4];

    if (size_t(start1)<start2) {
        pos[0] = size_t(start1);
        pos[1] = start2;
    }
    else {
        pos[0] = start2;
        pos[1] = size_t(start1);
    }

    pos[2] = helix->opposite_position(pos[1]);
    pos[3] = helix->opposite_position(pos[0]);

    return pos;
}

// -------------------------------------------
// current way to save folded helices:
//
// save a list of helix-numbers (with sign) separated by ';'
// strand start-positions are saved as 'num'
// segment start-positions are saved as '!num'

char *SEC_xstring_to_foldedHelixList(const char *x_string, size_t xlength, const BI_helix *helix, GB_ERROR& error) {
    GBS_strstruct out(1000);
    bool next_is_start = true;

    sec_assert(!error);
    error = NULp;

    for (size_t pos = 0; pos<xlength && !error; ++pos) {
        if (x_string[pos] == 'x') {
            const char *helixNr = helix->helixNr(pos);

            if (next_is_start) {
                sec_assert(helix->is_pairpos(pos));
                out.cat(helixNr);
            }
            else {
                bool is_end_and_start = helix->is_pairpos(pos);

                if (!helix->is_pairpos(pos-1)) {
                    // hackaround: xstring got out of sync (e.g. by insert/delete columns)
                    error = GBS_global_string("xstring out of sync at position %zu (approximately)\nRefold helices around position!", pos);
                }
                else {
                    out.nprintf(20, "!%s", helix->helixNr(pos-1));

                    if (is_end_and_start) {
                        out.put(';');
                        out.cat(helixNr);
                        next_is_start = !next_is_start;
                    }
                }
            }
            next_is_start = !next_is_start;
            out.put(';');
        }
    }

    return error ? NULp : out.release();
}

char *SEC_foldedHelixList_to_xstring(const char *foldedHelices, size_t xlength, const BI_helix *helix, GB_ERROR& error) {
    // example of foldedHelices: '-5;!-5;-8;!-8;-9;!-9;9;!9;8;!8;5;!5;-18;!-18;18;!18;-19a;!-19a;19a;!19a;-21;!-21;-24;!-24;24;!24;21;!21;-27;!-27;27;!27;-31;!-31;-43;!-43;43;!43;-46;!-46;46;!46;31;!31;-50;!-50;50;!50;'
    char *xstring    = ARB_alloc<char>(xlength+1);
    memset(xstring, '.', xlength);
    xstring[xlength] = 0;

    sec_assert(!error);
    error = NULp;

    while (!error) {
        const char *semi = strchr(foldedHelices, ';');
        if (!semi) break;

        char       *taggedHelixNr = ARB_strpartdup(foldedHelices, semi-1);
        const char *helixNr       = NULp;
        long        xpos          = -1;

        if (taggedHelixNr[0] == '!') { // position behind end of helix
            helixNr = taggedHelixNr+1;
            xpos    = helix->last_position(helixNr);
            if (xpos >= 0) ++xpos;
        }
        else { // position at start of helix
            helixNr = taggedHelixNr;
            xpos    = helix->first_position(helixNr);
        }

        if (xpos == -1) {
            error = GBS_global_string("Can't find helix '%s'", helixNr);
        }
        else {
            sec_assert(xpos >= 0 && size_t(xpos) < xlength);
            xstring[xpos] = 'x';
        }

        free(taggedHelixNr);

        foldedHelices = semi+1;
    }

    return xstring;
}

// -------------------------------------------------------------------------
//      xstring was made helix-relative, when saved to old version file


char *old_decode_xstring_rel_helix(GB_CSTR rel_helix, size_t xlength, const BI_helix *helix, int *no_of_helices_ptr) {
    const char *start_helix_nr = NULp;
    int         no_of_helices  = 0;
    int         start_helix    = 0;
    int         end_helix      = -1;
    int         rel_pos        = 0;
    size_t      lastpos        = helix->size()-1;

    char *x_buffer    = ARB_alloc<char>(xlength+1);
    memset(x_buffer, '.', xlength);
    x_buffer[xlength] = 0;

    for (size_t pos=0; ; pos++) {
        const char *helix_nr = NULp; // [required due to goto below]

        if (helix->is_pairpos(pos)) {
            helix_nr = helix->helixNr(pos);

            if (helix_nr==start_helix_nr) { // same helix as last
                end_helix = pos;
            }
            else { // new helix nr
                if (start_helix_nr) { // not first helix -> write last to xstring
                insert_helix :
                    helix_nr = helix->helixNr(pos); // re-init (needed in case of goto insert_helix)
                    char flag = rel_helix[rel_pos++];

                    no_of_helices++;
                    if (flag=='1') {
                        sec_assert(end_helix!=-1);
                        sec_assert(size_t(start_helix)<xlength);
                        sec_assert(size_t(end_helix+1)<xlength);

                        x_buffer[start_helix] = 'x';
                        x_buffer[end_helix+1] = 'x';
                    }
                    else if (flag!='0') {
                        if (flag==0) break; // eos

                        sec_assert(0); // illegal character

                        break;
                    }
                    if (pos==lastpos) {
                        break;
                    }
                }
                start_helix = pos;
                end_helix = pos;
                start_helix_nr = helix_nr;
            }
        }
        if (pos==lastpos) {
            if (start_helix_nr) {
                goto insert_helix;
            }
        }
    }

    *no_of_helices_ptr = no_of_helices;
    return x_buffer;
}