sequence.cpp

来自「ncbi源码」· C++ 代码 · 共 2,155 行 · 第 1/5 页

    // Check for overlap of mec and youc
    if ( got_one ) {
        return eOverlap;
    }
    mit = mec.begin();
    cmp1 = s_Compare(**mit, you, scope);
    for (++mit;  mit != mec.end();  ++mit) {
        cmp2 = s_Compare(**mit, you, scope);
        cmp1 = static_cast<ECompare> (s_RetA[cmp1][cmp2]);
    }
    return cmp1;
}


inline
const CSeq_loc::TLocations& s_GetLocations(const CSeq_loc& loc)
{
    switch (loc.Which()) {
        case CSeq_loc::e_Mix:    return loc.GetMix().Get();
        case CSeq_loc::e_Equiv:  return loc.GetEquiv().Get();
        default: // should never happen, but the compiler doesn't know that...
        NCBI_THROW(CObjmgrUtilException, eBadLocation,
                         "s_GetLocations: unsupported location type:"
                         + CSeq_loc::SelectionName(loc.Which()));
    }
}


// Compares two Seq-locs
ECompare s_Compare
(const CSeq_loc& me,
 const CSeq_loc& you,
 CScope*         scope)
{
    // Handle the case where me is one of e_Mix, e_Equiv, e_Packed_int
    switch ( me.Which() ) {
    case CSeq_loc::e_Mix:
    case CSeq_loc::e_Equiv: {
        const CSeq_loc::TLocations& pmc = s_GetLocations(me);
        switch ( you.Which() ) {
        case CSeq_loc::e_Mix:
        case CSeq_loc::e_Equiv: {
            const CSeq_loc::TLocations& pyc = s_GetLocations(you);
            return s_Compare_Help(pmc, pyc, you, scope);
        }
        case CSeq_loc::e_Packed_int: {
            const CSeq_loc::TIntervals& pyc = you.GetPacked_int().Get();
            return s_Compare_Help(pmc,pyc, you, scope);
        }
        case CSeq_loc::e_Null:
        case CSeq_loc::e_Empty:
        case CSeq_loc::e_Whole:
        case CSeq_loc::e_Int:
        case CSeq_loc::e_Pnt:
        case CSeq_loc::e_Packed_pnt:
        case CSeq_loc::e_Bond:
        case CSeq_loc::e_Feat: {
            //Check that pmc is not empty
            if(pmc.size() == 0) {
                return eNoOverlap;
            }

            CSeq_loc::TLocations::const_iterator mit = pmc.begin();
            ECompare cmp1, cmp2;
            cmp1 = s_Compare(**mit, you, scope);

            // This loop will only be executed if pmc.size() > 1
            for (++mit;  mit != pmc.end();  ++mit) {
                cmp2 = s_Compare(**mit, you, scope);
                cmp1 = static_cast<ECompare> (s_RetA[cmp1][cmp2]);
            }
            return cmp1;
        }
        default:
            return eNoOverlap;
        }
    }
    case CSeq_loc::e_Packed_int: {
        switch ( you.Which() ) {
        case CSeq_loc::e_Mix:
        case CSeq_loc::e_Equiv: {
            const CSeq_loc::TLocations& pyc = s_GetLocations(you);
            const CSeq_loc::TIntervals& pmc = me.GetPacked_int().Get();
            return s_Compare_Help(pmc,pyc, you, scope);
        }
        case CSeq_loc::e_Packed_int: {
            const CSeq_loc::TIntervals& mc = me.GetPacked_int().Get();
            const CSeq_loc::TIntervals& yc = you.GetPacked_int().Get();
            return s_Compare_Help(mc, yc, you, scope);
        }
        default:
            // All other cases are handled below
            break;
        }
    }
    default:
        // All other cases are handled below
        break;
    }

    // Note, at this point, me is not one of e_Mix or e_Equiv
    switch ( you.Which() ) {
    case CSeq_loc::e_Mix:
    case CSeq_loc::e_Equiv: {
        const CSeq_loc::TLocations& pyouc = s_GetLocations(you);

        // Check that pyouc is not empty
        if(pyouc.size() == 0) {
            return eNoOverlap;
        }

        CSeq_loc::TLocations::const_iterator yit = pyouc.begin();
        ECompare cmp1, cmp2;
        cmp1 = s_Compare(me, **yit, scope);

        // This loop will only be executed if pyouc.size() > 1
        for (++yit;  yit != pyouc.end();  ++yit) {
            cmp2 = s_Compare(me, **yit, scope);
            cmp1 = static_cast<ECompare> (s_RetB[cmp1][cmp2]);
        }
        return cmp1;
    }
    // All other cases handled below
    default:
        break;
    }

    switch ( me.Which() ) {
    case CSeq_loc::e_Null: {
        switch ( you.Which() ) {
        case CSeq_loc::e_Null:
            return eSame;
        default:
            return eNoOverlap;
        }
    }
    case CSeq_loc::e_Empty: {
        switch ( you.Which() ) {
        case CSeq_loc::e_Empty:
            if ( IsSameBioseq(me.GetEmpty(), you.GetEmpty(), scope) ) {
                return eSame;
            } else {
                return eNoOverlap;
            }
        default:
            return eNoOverlap;
        }
    }
    case CSeq_loc::e_Packed_int: {
        // Comparison of two e_Packed_ints is handled above
        switch ( you.Which() ) {
        case CSeq_loc::e_Whole:
            return s_Compare(me.GetPacked_int(), you.GetWhole(), scope);
        case CSeq_loc::e_Int:
            return s_Compare(me.GetPacked_int(), you.GetInt(), scope);
        case CSeq_loc::e_Pnt:
            return s_Compare(me.GetPacked_int(), you.GetPnt(), scope);
        case CSeq_loc::e_Packed_pnt:
            return s_Compare(me.GetPacked_int(), you.GetPacked_pnt(), scope);
        case CSeq_loc::e_Bond:
            return s_Compare(me.GetPacked_int(), you.GetBond(), scope);
        default:
            return eNoOverlap;
        }
    }
    case CSeq_loc::e_Whole: {
        switch ( you.Which() ) {
        case CSeq_loc::e_Whole:
            return s_Compare(me.GetWhole(), you.GetWhole(), scope);
        case CSeq_loc::e_Bond:
            return s_Compare(me.GetWhole(), you.GetBond(), scope);
        case CSeq_loc::e_Int:
            return s_Compare(me.GetWhole(), you.GetInt(), scope);
        case CSeq_loc::e_Pnt:
            return s_Compare(me.GetWhole(), you.GetPnt(), scope);
        case CSeq_loc::e_Packed_pnt:
            return s_Compare(me.GetWhole(), you.GetPacked_pnt(), scope);
        case CSeq_loc::e_Packed_int:
            return s_Compare(me.GetWhole(), you.GetPacked_int(), scope);
        default:
            return eNoOverlap;
        }
    }
    case CSeq_loc::e_Int: {
        switch ( you.Which() ) {
        case CSeq_loc::e_Whole:
            return s_Compare(me.GetInt(), you.GetWhole(), scope);
        case CSeq_loc::e_Bond:
            return s_Compare(me.GetInt(), you.GetBond(), scope);
        case CSeq_loc::e_Int:
            return s_Compare(me.GetInt(), you.GetInt(), scope);
        case CSeq_loc::e_Pnt:
            return s_Compare(me.GetInt(), you.GetPnt(), scope);
        case CSeq_loc::e_Packed_pnt:
            return s_Compare(me.GetInt(), you.GetPacked_pnt(), scope);
        case CSeq_loc::e_Packed_int:
            return s_Compare(me.GetInt(), you.GetPacked_int(), scope);
        default:
            return eNoOverlap;
        }
    }
    case CSeq_loc::e_Pnt: {
        switch ( you.Which() ) {
        case CSeq_loc::e_Whole:
            return s_Complement(s_Compare(you.GetWhole(), me.GetPnt(), scope));
        case CSeq_loc::e_Int:
            return s_Compare(me.GetPnt(), you.GetInt(), scope);
        case CSeq_loc::e_Pnt:
            return s_Compare(me.GetPnt(), you.GetPnt(), scope);
        case CSeq_loc::e_Packed_pnt:
            return s_Compare(me.GetPnt(), you.GetPacked_pnt(), scope);
        case CSeq_loc::e_Bond:
            return s_Compare(me.GetPnt(), you.GetBond(), scope);
        case CSeq_loc::e_Packed_int:
            return s_Compare(me.GetPnt(), you.GetPacked_int(), scope);
        default:
            return eNoOverlap;
        }
    }
    case CSeq_loc::e_Packed_pnt: {
        const CPacked_seqpnt& packed = me.GetPacked_pnt();
        switch ( you.Which() ) {
        case CSeq_loc::e_Whole:
            return s_Complement(s_Compare(you.GetWhole(), packed, scope));
        case CSeq_loc::e_Int:
            return s_Compare(packed, you.GetInt(), scope);
        case CSeq_loc::e_Pnt:
            return s_Complement(s_Compare(you.GetPnt(), packed, scope));
        case CSeq_loc::e_Packed_pnt:
            return s_Compare(packed, you.GetPacked_pnt(),scope);
        case CSeq_loc::e_Bond:
            return s_Compare(packed, you.GetBond(), scope);
        case CSeq_loc::e_Packed_int:
            return s_Compare(me.GetPacked_pnt(), you.GetPacked_int(), scope);
        default:
            return eNoOverlap;
        }
    }
    case CSeq_loc::e_Bond: {
        const CSeq_bond& bnd = me.GetBond();
        switch ( you.Which() ) {
        case CSeq_loc::e_Whole:
            return s_Complement(s_Compare(you.GetWhole(), bnd, scope));
        case CSeq_loc::e_Bond:
            return s_Compare(bnd, you.GetBond(), scope);
        case CSeq_loc::e_Int:
            return s_Compare(bnd, you.GetInt(), scope);
        case CSeq_loc::e_Packed_pnt:
            return s_Complement(s_Compare(you.GetPacked_pnt(), bnd, scope));
        case CSeq_loc::e_Pnt:
            return s_Complement(s_Compare(you.GetPnt(), bnd, scope));
        case CSeq_loc::e_Packed_int:
            return s_Complement(s_Compare(you.GetPacked_int(), bnd, scope));
        default:
            return eNoOverlap;
        }
    }
    default:
        return eNoOverlap;
    }
}


ECompare Compare
(const CSeq_loc& loc1,
 const CSeq_loc& loc2,
 CScope*         scope)
{
    return s_Compare(loc1, loc2, scope);
}


void ChangeSeqId(CSeq_id* id, bool best, CScope* scope)
{
    // Return if no scope
    if (!scope) {
        return;
    }

    // Get CBioseq represented by *id
    CBioseq_Handle::TBioseqCore seq =
        scope->GetBioseqHandle(*id).GetBioseqCore();

    // Get pointer to the best/worst id of *seq
    const CSeq_id* tmp_id;
    if (!best) {
        tmp_id = FindBestChoice(seq->GetId(), CSeq_id::BestRank).GetPointer();
    } else {
        tmp_id = FindBestChoice(seq->GetId(), CSeq_id::WorstRank).GetPointer();
    }

    // Change the contents of *id to that of *tmp_id
    id->Reset();
    id->Assign(*tmp_id);
}


void ChangeSeqLocId(CSeq_loc* loc, bool best, CScope* scope)
{
    if (!scope) {
        return;
    }

    for (CTypeIterator<CSeq_id> id(Begin(*loc)); id; ++id) {
        ChangeSeqId(&(*id), best, scope);
    }
}


bool BadSeqLocSortOrder
(const CBioseq&, //  seq,
 const CSeq_loc& loc,
 CScope*         scope)
{
    ENa_strand strand = GetStrand(loc, scope);
    if (strand == eNa_strand_unknown  ||  strand == eNa_strand_other) {
        return false;
    }
    
    // Check that loc segments are in order
    CSeq_loc::TRange last_range;
    bool first = true;
    for (CSeq_loc_CI lit(loc); lit; ++lit) {
        if (first) {
            last_range = lit.GetRange();
            first = false;
            continue;
        }
        if (strand == eNa_strand_minus) {
            if (last_range.GetTo() < lit.GetRange().GetTo()) {
                return true;
            }
        } else {
            if (last_range.GetFrom() > lit.GetRange().GetFrom()) {
                return true;
            }
        }
        last_range = lit.GetRange();
    }
    return false;
}


ESeqLocCheck SeqLocCheck(const CSeq_loc& loc, CScope* scope)
{
    ESeqLocCheck rtn = eSeqLocCheck_ok;

    ENa_strand strand = GetStrand(loc, scope);
    if (strand == eNa_strand_unknown  ||  strand == eNa_strand_other) {
        rtn = eSeqLocCheck_warning;
    }

    CTypeConstIterator<CSeq_loc> lit(ConstBegin(loc));
    for (;lit; ++lit) {
        switch (lit->Which()) {
        case CSeq_loc::e_Int:
            if (!IsValid(lit->GetInt(), scope)) {
                rtn = eSeqLocCheck_error;
            }
            break;
        case CSeq_loc::e_Packed_int:
        {
            CTypeConstIterator<CSeq_interval> sit(ConstBegin(*lit));
            for(;sit; ++sit) {
                if (!IsValid(*sit, scope)) {
                    rtn = eSeqLocCheck_error;
                    break;
                }
            }
            break;
        }
        case CSeq_loc::e_Pnt:
            if (!IsValid(lit->GetPnt(), scope)) {
                rtn = eSeqLocCheck_error;
            }
            break;
        case CSeq_loc::e_Packed_pnt:
            if (!IsValid(lit->GetPacked_pnt(), scope)) {
                rtn = eSeqLocCheck_error;
            }
            break;
        default:
            break;
        }
    }
    return rtn;
}


CRef<CSeq_loc> SourceToProduct(const CSeq_feat& feat,
                               const CSeq_loc& source_loc, TS2PFlags flags,
                               CScope* scope, int* frame)
{
    SRelLoc::TFlags rl_flags = 0;
    if (flags & fS2P_NoMerge) {
        rl_flags |= SRelLoc::fNoMerge;
    }
    SRelLoc rl(feat.GetLocation(), source_loc, scope, rl_flags);
    _ASSERT(!rl.m_Ranges.empty());
    rl.m_ParentLoc.Reset(&feat.GetProduct());
    if (feat.GetData().IsCdregion()) {
        // 3:1 ratio
        const CCdregion& cds         = feat.GetData().GetCdregion();
        int              base_frame  = cds.GetFrame();
        if (base_frame > 0) {
            --base_frame;
        }
        if (frame) {
            *frame = 3 - (rl.m_Ranges.front()->GetFrom() + 2 - base_frame) % 3;
        }
        TSeqPos prot_length;
        try {
            prot_length = GetLength(feat.GetProduct(), scope);
        } catch (CNoLength) {
            prot_length = numeric_limits<TSeqPos>::max();
        }
        NON_CONST_ITERATE (SRelLoc::TRanges, it, rl.m_Ranges) {
            if (IsReverse((*it)->GetStrand())) {
                ERR_POST(Warning
                         << "SourceToProduct:"
                         " parent and child have opposite orientations");
            }
            (*it)->SetFrom(((*it)->GetFrom() - base_frame) / 3);
            (*it)->SetTo  (((*it)->GetTo()   - base_frame) / 3);
            if ((flags & fS2P_AllowTer)  &&  (*it)->GetTo() == prot_length) {
                --(*it)->SetTo();
            }
        }
    } else {
        if (frame) {
            *frame = 0; // not applicable; explicitly zero