sequence.cpp

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

                strand_set = true;
            } else if (!strand_set) {
                strand = istrand;
                strand_set = true;
            } else if (istrand != strand) {
                return eNa_strand_other;
            }
        }
    }
    default:
        return eNa_strand_unknown;
    }
}


ENa_strand GetStrand(const CSeq_loc& loc, CScope* scope)
{
    if (!IsOneBioseq(loc, scope)) {
        return eNa_strand_unknown;  // multiple bioseqs
    }

    ENa_strand strand = eNa_strand_unknown, cstrand;
    bool strand_set = false;
    for (CSeq_loc_CI i(loc); i; ++i) {
        switch (i.GetSeq_loc().Which()) {
        case CSeq_loc::e_Equiv:
            break;
        default:
            cstrand = s_GetStrand(i.GetSeq_loc());
            if (strand == eNa_strand_unknown  &&  cstrand == eNa_strand_plus) {
                strand = eNa_strand_plus;
                strand_set = true;
            } else if (strand == eNa_strand_plus  &&
                cstrand == eNa_strand_unknown) {
                cstrand = eNa_strand_plus;
                strand_set = true;
            } else if (!strand_set) {
                strand = cstrand;
                strand_set = true;
            } else if (cstrand != strand) {
                return eNa_strand_other;
            }
        }
    }
    return strand;
}


TSeqPos GetStart(const CSeq_loc& loc, CScope* scope)
    THROWS((CNotUnique))
{
    // Throw CNotUnique if loc does not represent one CBioseq
    GetId(loc, scope);

    CSeq_loc::TRange rg = loc.GetTotalRange();
    return rg.GetFrom();
}


TSeqPos GetStop(const CSeq_loc& loc, CScope* scope)
    THROWS((CNotUnique))
{
    // Throw CNotUnique if loc does not represent one CBioseq
    GetId(loc, scope);

    CSeq_loc::TRange rg = loc.GetTotalRange();
    return rg.GetTo();
}


/*
*******************************************************************************
                        Implementation of Compare()
*******************************************************************************
*/


ECompare s_Compare
(const CSeq_loc&,
 const CSeq_loc&,
 CScope*);

int s_RetA[5][5] = {
    { 0, 4, 2, 2, 4 },
    { 4, 1, 4, 1, 4 },
    { 2, 4, 2, 2, 4 },
    { 2, 1, 2, 3, 4 },
    { 4, 4, 4, 4, 4 }
};


int s_RetB[5][5] = {
    { 0, 1, 4, 1, 4 },
    { 1, 1, 4, 1, 1 },
    { 4, 4, 2, 2, 4 },
    { 1, 1, 4, 3, 4 },
    { 4, 1, 4, 4, 4 }
};


// Returns the complement of an ECompare value
inline
ECompare s_Complement(ECompare cmp)
{
    switch ( cmp ) {
    case eContains:
        return eContained;
    case eContained:
        return eContains;
    default:
        return cmp;
    }
}


// Compare two Whole sequences
inline
ECompare s_Compare
(const CSeq_id& id1,
 const CSeq_id& id2,
 CScope*        scope)
{
    // If both sequences from same CBioseq, they are the same
    if ( IsSameBioseq(id1, id2, scope) ) {
        return eSame;
    } else {
        return eNoOverlap;
    }
}


// Compare Whole sequence with a Bond
inline
ECompare s_Compare
(const CSeq_id&   id,
 const CSeq_bond& bond,
 CScope*          scope)
{
    unsigned int count = 0;

    // Increment count if bond point A is from same CBioseq as id
    if ( IsSameBioseq(id, bond.GetA().GetId(), scope) ) {
        ++count;
    }

    // Increment count if second point in bond is set and is from
    // same CBioseq as id.
    if (bond.IsSetB()  &&  IsSameBioseq(id, bond.GetB().GetId(), scope)) {
        ++count;
    }

    switch ( count ) {
    case 1:
        if ( bond.IsSetB() ) {
            // One of two bond points are from same CBioseq as id --> overlap
            return eOverlap;
        } else {
            // There is only one bond point set --> id contains bond
            return eContains;
        }
    case 2:
        // Both bond points are from same CBioseq as id --> id contains bond
        return eContains;
    default:
        // id and bond do not overlap
        return eNoOverlap;
    }
}


// Compare Whole sequence with an interval
inline
ECompare s_Compare
(const CSeq_id&       id,
 const CSeq_interval& interval,
 CScope*              scope)
{
    // If interval is from same CBioseq as id and interval is
    // [0, length of id-1], then they are the same. If interval is from same
    // CBioseq as id, but interval is not [0, length of id -1] then id
    // contains seqloc.
    if ( IsSameBioseq(id, interval.GetId(), scope) ) {
        if (interval.GetFrom() == 0  &&
            interval.GetTo()   == GetLength(id, scope) - 1) {
            return eSame;
        } else {
            return eContains;
        }
    } else {
        return eNoOverlap;
    }
}


// Compare Whole sequence with a point
inline
ECompare s_Compare
(const CSeq_id&     id,
 const CSeq_point&  point,
 CScope*            scope)
{
    // If point from the same CBioseq as id, then id contains point
    if ( IsSameBioseq(id, point.GetId(), scope) ) {
        return eContains;
    } else {
        return eNoOverlap;
    }
}


// Compare Whole sequence with packed points
inline
ECompare s_Compare
(const CSeq_id&        id,
 const CPacked_seqpnt& packed,
 CScope*               scope)
{
    // If packed from same CBioseq as id, then id contains packed
    if ( IsSameBioseq(id, packed.GetId(), scope) ) {
        return eContains;
    } else {
        return eNoOverlap;
    }
}


// Compare an interval with a bond
inline
ECompare s_Compare
(const CSeq_interval& interval,
 const CSeq_bond&     bond,
 CScope*              scope)
{
    unsigned int count = 0;

    // Increment count if interval contains the first point of the bond
    if (IsSameBioseq(interval.GetId(), bond.GetA().GetId(), scope)  &&
        interval.GetFrom() <= bond.GetA().GetPoint()  &&
        interval.GetTo()   >= bond.GetA().GetPoint()) {
        ++count;
    }

    // Increment count if the second point of the bond is set and the
    // interval contains the second point.
    if (bond.IsSetB()  &&
        IsSameBioseq(interval.GetId(), bond.GetB().GetId(), scope)  &&
        interval.GetFrom() <= bond.GetB().GetPoint()  &&
        interval.GetTo()   >= bond.GetB().GetPoint()) {
        ++count;
    }

    switch ( count ) {
    case 1:
        if ( bond.IsSetB() ) {
            // The 2nd point of the bond is set
            return eOverlap;
        } else {
            // The 2nd point of the bond is not set
            return eContains;
        }
    case 2:
        // Both points of the bond are in the interval
        return eContains;
    default:
        return eNoOverlap;
    }
}


// Compare a bond with an interval
inline
ECompare s_Compare
(const CSeq_bond&     bond,
 const CSeq_interval& interval,
 CScope*              scope)
{
    // Just return the complement of comparing an interval with a bond
    return s_Complement(s_Compare(interval, bond, scope));
}


// Compare an interval to an interval
inline
ECompare s_Compare
(const CSeq_interval& intA,
 const CSeq_interval& intB,
 CScope*              scope)
{
    // If the intervals are not on the same bioseq, then there is no overlap
    if ( !IsSameBioseq(intA.GetId(), intB.GetId(), scope) ) {
        return eNoOverlap;
    }

    // Compare two intervals
    TSeqPos fromA = intA.GetFrom();
    TSeqPos toA   = intA.GetTo();
    TSeqPos fromB = intB.GetFrom();
    TSeqPos toB   = intB.GetTo();

    if (fromA == fromB  &&  toA == toB) {
        // End points are the same --> the intervals are the same.
        return eSame;
    } else if (fromA <= fromB  &&  toA >= toB) {
        // intA contains intB
        return eContains;
    } else if (fromA >= fromB  &&  toA <= toB) {
        // intB contains intA
        return eContained;
    } else if (fromA > toB  ||  toA < fromB) {
        // No overlap
        return eNoOverlap;
    } else {
        // The only possibility left is overlap
        return eOverlap;
    }
}


// Compare an interval and a point
inline
ECompare s_Compare
(const CSeq_interval& interval,
 const CSeq_point&    point,
 CScope*              scope)
{
    // If not from same CBioseq, then no overlap
    if ( !IsSameBioseq(interval.GetId(), point.GetId(), scope) ) {
        return eNoOverlap;
    }

    TSeqPos pnt = point.GetPoint();

    // If the interval is of length 1 and contains the point, then they are
    // identical
    if (interval.GetFrom() == pnt  &&  interval.GetTo() == pnt ) {
        return eSame;
    }

    // If point in interval, then interval contains point
    if (interval.GetFrom() <= pnt  &&  interval.GetTo() >= pnt ) {
        return eContains;
    }

    // Only other possibility
    return eNoOverlap;
}


// Compare a point and an interval
inline
ECompare s_Compare
(const CSeq_point&    point,
 const CSeq_interval& interval,
 CScope*              scope)
{
    // The complement of comparing an interval with a point.
    return s_Complement(s_Compare(interval, point, scope));
}


// Compare a point with a point
inline
ECompare s_Compare
(const CSeq_point& pntA,
 const CSeq_point& pntB,
 CScope*           scope)
{
    // If points are on same bioseq and coordinates are the same, then same.
    if (IsSameBioseq(pntA.GetId(), pntB.GetId(), scope)  &&
        pntA.GetPoint() == pntB.GetPoint() ) {
        return eSame;
    }

    // Otherwise they don't overlap
    return eNoOverlap;
}


// Compare a point with packed point
inline
ECompare s_Compare
(const CSeq_point&     point,
 const CPacked_seqpnt& points,
 CScope*               scope)
{
    // If on the same bioseq, and any of the packed points are the
    // same as point, then the point is contained in the packed point
    if ( IsSameBioseq(point.GetId(), points.GetId(), scope) ) {
        TSeqPos pnt = point.GetPoint();

        // This loop will only be executed if points.GetPoints().size() > 0
        ITERATE(CSeq_loc::TPoints, it, points.GetPoints()) {
            if (pnt == *it) {
                return eContained;
            }
        }
    }

    // Otherwise, no overlap
    return eNoOverlap;
}


// Compare a point with a bond
inline
ECompare s_Compare
(const CSeq_point& point,
 const CSeq_bond&  bond,
 CScope*           scope)
{
    unsigned int count = 0;

    // If the first point of the bond is on the same CBioseq as point
    // and the point coordinates are the same, increment count by 1
    if (IsSameBioseq(point.GetId(), bond.GetA().GetId(), scope)  &&
        bond.GetA().GetPoint() == point.GetPoint()) {
        ++count;
    }

    // If the second point of the bond is set, the points are from the
    // same CBioseq, and the coordinates of the points are the same,
    // increment the count by 1.
    if (bond.IsSetB()  &&
        IsSameBioseq(point.GetId(), bond.GetB().GetId(), scope)  &&
        bond.GetB().GetPoint() == point.GetPoint()) {
        ++count;
    }

    switch ( count ) {
    case 1:
        if ( bond.IsSetB() ) {