nw_formatter.cpp

ncbi源码

                line3[jPos] = c2;
            }
            ss << line2 << endl << line3 << endl << endl;
        }
    }
    break;

    case eFormatAsn: {
        CSeq_align seq_align;
        AsSeqAlign(&seq_align);
        CObjectOStreamAsn asn_stream (ss);
        asn_stream << seq_align;
        asn_stream << Separator;
    }
    break;

    case eFormatFastA: {
        vector<char> v1, v2;
        size_t aln_size = x_ApplyTranscript(&v1, &v2);
        
        ss << '>' << m_Seq1Id << endl;
        const vector<char>* pv = &v1;
        for(size_t i = 0; i < aln_size; ++i) {
            for(size_t j = 0; j < line_width && i < aln_size; ++j, ++i) {
                ss << (*pv)[i];
            }
            ss << endl;
        }

        ss << '>' << m_Seq2Id << endl;
        pv = &v2;
        for(size_t i = 0; i < aln_size; ++i) {
            for(size_t j = 0; j < line_width && i < aln_size; ++j, ++i) {
                ss << (*pv)[i];
            }
            ss << endl;
        }
    }
    break;
    
    case eFormatExonTable:
    case eFormatExonTableEx: {
        ss.precision(3);
        bool esfL1, esfR1, esfL2, esfR2;
        m_aligner->GetEndSpaceFree(&esfL1, &esfR1, &esfL2, &esfR2);
        const size_t len2 = m_aligner->GetSeqLen2();
        const char* start1 = m_aligner->GetSeq1();
        const char* start2 = m_aligner->GetSeq2();
        const char* p1 = start1;
        const char* p2 = start2;
        int tr_idx_hi0 = transcript.size() - 1, tr_idx_hi = tr_idx_hi0;
        int tr_idx_lo0 = 0, tr_idx_lo = tr_idx_lo0;

        if(esfL1 && transcript[tr_idx_hi0] == CNWAligner::eTS_Insert) {
            while(esfL1 && transcript[tr_idx_hi] == CNWAligner::eTS_Insert) {
                --tr_idx_hi;
                ++p2;
            }
        }

        if(esfL2 && transcript[tr_idx_hi0] == CNWAligner::eTS_Delete) {
            while(esfL2 && transcript[tr_idx_hi] == CNWAligner::eTS_Delete) {
                --tr_idx_hi;
                ++p1;
            }
        }

        if(esfR1 && transcript[tr_idx_lo0] == CNWAligner::eTS_Insert) {
            while(esfR1 && transcript[tr_idx_lo] == CNWAligner::eTS_Insert) {
                ++tr_idx_lo;
            }
        }

        if(esfR2 && transcript[tr_idx_lo0] == CNWAligner::eTS_Delete) {
            while(esfR2 && transcript[tr_idx_lo] == CNWAligner::eTS_Delete) {
                ++tr_idx_lo;
            }
        }

        bool type_ex = type == eFormatExonTableEx;
        vector<char> trans_ex (tr_idx_hi - tr_idx_lo + 1);

        for(int tr_idx = tr_idx_hi; tr_idx >= tr_idx_lo;) {
            const char* p1_beg = p1;
            const char* p2_beg = p2;
            size_t matches = 0, exon_aln_size = 0;

            vector<char>::iterator ii_ex = trans_ex.begin();
            while(tr_idx >= tr_idx_lo && 
                  transcript[tr_idx] < CNWAligner::eTS_Intron) {
                
                bool noins = transcript[tr_idx] != CNWAligner::eTS_Insert;
                bool nodel = transcript[tr_idx] != CNWAligner::eTS_Delete;
                if(noins && nodel) {
                    if(*p1++ == *p2++) {
                        ++matches;
                        if(type_ex) *ii_ex++ = 'M';
                    }
                    else {
                        if(type_ex) *ii_ex++ = 'R';
                    }
                } else if( noins ) {
                    ++p1;
                    if(type_ex) *ii_ex++ = 'D';
                } else {
                    ++p2;
                    if(type_ex) *ii_ex++ = 'I';
                }
                --tr_idx;
                ++exon_aln_size;
            }

            if(exon_aln_size > 0) {

                if(m_Seq1Id.size() && m_Seq2Id.size()) {
                    ss << m_Seq1Id << '\t' << m_Seq2Id << '\t';
                }
		else {
		    ss << "-\t-\t";
		}
                float identity = float(matches) / exon_aln_size;
                ss << identity << '\t' << exon_aln_size << '\t';
                size_t beg1  = p1_beg - start1, end1 = p1 - start1 - 1;
                size_t beg2  = p2_beg - start2, end2 = p2 - start2 - 1;
		if(beg1 <= end1) {
		    ss << beg1 << '\t' << end1 << '\t';
		}
		else {
		    ss << "-\t-\t";
		}
		ss << beg2 << '\t' << end2 << '\t';
                char c1 = (p2_beg >= start2 + 2)? *(p2_beg - 2): ' ';
                char c2 = (p2_beg >= start2 + 1)? *(p2_beg - 1): ' ';
                char c3 = (p2 < start2 + len2)? *(p2): ' ';
                char c4 = (p2 < start2 + len2 - 1)? *(p2+1): ' ';
                ss << c1 << c2 << "<exon>" << c3 << c4;

                if(type == eFormatExonTableEx) {
                    ss << '\t';
                    copy(trans_ex.begin(), ii_ex, ostream_iterator<char>(ss));
                }
                ss << endl;
            }

            // find next exon
            while(tr_idx >= tr_idx_lo &&
                  (transcript[tr_idx] == CNWAligner::eTS_Intron)) {

                --tr_idx;
                ++p2;
            }
        }
    }
    break;

    default: {
        NCBI_THROW(CAlgoAlignException, eBadParameter,
                   "Incorrect format specified");
    }

    }

    *output = CNcbiOstrstreamToString(ss);
}



// Transform source sequences according to the transcript.
// Write the results to v1 and v2 leaving source sequences intact.
// Return alignment size.
size_t CNWFormatter::x_ApplyTranscript(vector<char>* pv1, vector<char>* pv2)
    const
{
    const vector<CNWAligner::ETranscriptSymbol>& transcript = 
        *(m_aligner->GetTranscript());

    vector<char>& v1 = *pv1;
    vector<char>& v2 = *pv2;

    vector<CNWAligner::ETranscriptSymbol>::const_reverse_iterator
        ib = transcript.rbegin(),
        ie = transcript.rend(),
        ii;

    const char* iv1 = m_aligner->GetSeq1();
    const char* iv2 = m_aligner->GetSeq2();
    v1.clear();
    v2.clear();

    for (ii = ib;  ii != ie;  ii++) {
        CNWAligner::ETranscriptSymbol ts = *ii;
        char c1, c2;
        switch ( ts ) {
        case CNWAligner::eTS_Insert:
            c1 = '-';
            c2 = *iv2++;
            break;
        case CNWAligner::eTS_Delete:
            c2 = '-';
            c1 = *iv1++;
            break;
        case CNWAligner::eTS_Match:
        case CNWAligner::eTS_Replace:
            c1 = *iv1++;
            c2 = *iv2++;
            break;
        case CNWAligner::eTS_Intron:
            c1 = '+';
            c2 = *iv2++;
            break;
        default:
            c1 = c2 = '?';
            break;
        }
        v1.push_back(c1);
        v2.push_back(c2);
    }

    return v1.size();
}


END_NCBI_SCOPE


/*
 * ===========================================================================
 * $Log: nw_formatter.cpp,v $
 * Revision 1000.2  2004/06/01 18:04:54  gouriano
 * PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.10
 *
 * Revision 1.10  2004/05/21 21:41:02  gorelenk
 * Added PCH ncbi_pch.hpp
 *
 * Revision 1.9  2004/05/18 21:43:40  kapustin
 * Code cleanup
 *
 * Revision 1.8  2004/05/17 14:50:56  kapustin
 * Add/remove/rearrange some includes and object declarations
 *
 * Revision 1.7  2004/03/18 16:30:24  grichenk
 * Changed type of seq-align containers from list to vector.
 *
 * Revision 1.6  2003/10/14 19:26:59  kapustin
 * Format void exons properly
 *
 * Revision 1.5  2003/09/26 14:43:18  kapustin
 * Remove exception specifications
 *
 * Revision 1.4  2003/09/15 21:32:03  kapustin
 * Minor code cleanup
 *
 * Revision 1.3  2003/09/12 19:43:04  kapustin
 * Add checking for empty transcript
 *
 * Revision 1.2  2003/09/03 01:19:32  ucko
 * +<iterator> (needed for ostream_iterator<> with some compilers)
 *
 * Revision 1.1  2003/09/02 22:34:49  kapustin
 * Initial revision
 *
 * ===========================================================================
 */