su_block_multiple_alignment.cpp

ncbi源码

        alignedBlock = *a;
        newUnalignedBlock = CreateNewUnalignedBlockBetween(prevAlignedBlock, alignedBlock);
        if (newUnalignedBlock)
            m_blocks.insert(a, CRef<Block>(newUnalignedBlock));
        prevAlignedBlock = alignedBlock;
    }

    // right tail
    newUnalignedBlock = CreateNewUnalignedBlockBetween(alignedBlock, NULL);
    if (newUnalignedBlock) {
        m_blocks.insert(a, CRef<Block>(newUnalignedBlock));
    }

    return true;
}

bool BlockMultipleAlignment::UpdateBlockMap(bool clearRowInfo)
{
    unsigned int i = 0, j, n = 0;
    BlockList::iterator b, be = m_blocks.end();

    // reset old stuff, recalculate m_width
    m_totalWidth = 0;
    for (b=m_blocks.begin(); b!=be; ++b)
        m_totalWidth += (*b)->m_width;

    // fill out the block map
    m_blockMap.resize(m_totalWidth);
    UngappedAlignedBlock *aBlock;
    for (b=m_blocks.begin(); b!=be; ++b) {
        aBlock = dynamic_cast<UngappedAlignedBlock*>(b->GetPointer());
        if (aBlock)
            ++n;
        for (j=0; j<(*b)->m_width; ++j, ++i) {
            m_blockMap[i].block = *b;
            m_blockMap[i].blockColumn = j;
            m_blockMap[i].alignedBlockNum = aBlock ? n : eUndefined;
        }
    }

    // if alignment changes, any pssm/scores/status become invalid
    RemovePSSM();
    if (clearRowInfo) ClearRowInfo();

    return true;
}

bool BlockMultipleAlignment::GetCharacterAt(
    unsigned int alignmentColumn, unsigned int row, eUnalignedJustification justification,
    char *character) const
{
    const Sequence *sequence;
    unsigned int seqIndex;
    bool isAligned;

    if (!GetSequenceAndIndexAt(alignmentColumn, row, justification, &sequence, &seqIndex, &isAligned))
        return false;

    *character = (seqIndex >= 0) ? sequence->m_sequenceString[seqIndex] : '~';
    if (isAligned)
        *character = toupper(*character);
    else
        *character = tolower(*character);

    return true;
}

bool BlockMultipleAlignment::GetSequenceAndIndexAt(
    unsigned int alignmentColumn, unsigned int row, eUnalignedJustification requestedJustification,
    const Sequence **sequence, unsigned int *index, bool *isAligned) const
{
    if (sequence)
        *sequence = m_sequences[row];

    const BlockInfo& blockInfo = m_blockMap[alignmentColumn];

    if (!blockInfo.block->IsAligned()) {
        if (isAligned)
            *isAligned = false;
        // override requested justification for end m_blocks
        if (blockInfo.block == m_blocks.back().GetPointer()) // also true if there's a single aligned block
            requestedJustification = eLeft;
        else if (blockInfo.block == m_blocks.front().GetPointer())
            requestedJustification = eRight;
    } else
        if (isAligned)
            *isAligned = true;

    if (index)
        *index = blockInfo.block->GetIndexAt(blockInfo.blockColumn, row, requestedJustification);

    return true;
}

bool BlockMultipleAlignment::IsAligned(unsigned int row, unsigned int seqIndex) const
{
    const Block *block = GetBlock(row, seqIndex);
    return (block && block->IsAligned());
}

const Block * BlockMultipleAlignment::GetBlock(unsigned int row, unsigned int seqIndex) const
{
    // make sure we're in range for this sequence
    if (row < 0 || seqIndex < 0 || row >= NRows() || seqIndex >= m_sequences[row]->Length()) {
        ERROR_MESSAGE("BlockMultipleAlignment::GetBlock() - coordinate out of range");
        return NULL;
    }

    const Block::Range *range;

    // first check to see if it's in the same block as last time.
    if (m_cachePrevBlock) {
        range = m_cachePrevBlock->GetRangeOfRow(row);
        if (seqIndex >= range->from && seqIndex <= range->to)
            return m_cachePrevBlock;
        ++m_cacheBlockIterator; // start search at next block
    } else {
        m_cacheBlockIterator = m_blocks.begin();
    }

    // otherwise, perform block search. This search is most efficient when queries
    // happen in order from left to right along a given row.
    do {
        if (m_cacheBlockIterator == m_blocks.end())
            m_cacheBlockIterator = m_blocks.begin();
        range = (*m_cacheBlockIterator)->GetRangeOfRow(row);
        if (seqIndex >= range->from && seqIndex <= range->to) {
            m_cachePrevBlock = *m_cacheBlockIterator; // cache this block
            return m_cachePrevBlock;
        }
        ++m_cacheBlockIterator;
    } while (1);
}

unsigned int BlockMultipleAlignment::GetFirstAlignedBlockPosition(void) const
{
    BlockList::const_iterator b = m_blocks.begin();
    if (m_blocks.size() > 0 && (*b)->IsAligned()) // first block is aligned
        return 0;
    else if (m_blocks.size() >= 2 && (*(++b))->IsAligned()) // second block is aligned
        return m_blocks.front()->m_width;
    else
        return eUndefined;
}

unsigned int BlockMultipleAlignment::GetAlignedSlaveIndex(unsigned int masterSeqIndex, unsigned int slaveRow) const
{
    const UngappedAlignedBlock
        *aBlock = dynamic_cast<const UngappedAlignedBlock*>(GetBlock(0, masterSeqIndex));
    if (!aBlock)
        return eUndefined;

    const Block::Range
        *masterRange = aBlock->GetRangeOfRow(0),
        *slaveRange = aBlock->GetRangeOfRow(slaveRow);
    return (slaveRange->from + masterSeqIndex - masterRange->from);
}

void BlockMultipleAlignment::GetAlignedBlockPosition(unsigned int alignmentIndex,
    unsigned int *blockColumn, unsigned int *blockWidth) const
{
    *blockColumn = *blockWidth = eUndefined;
    const BlockInfo& info = m_blockMap[alignmentIndex];
    if (info.block->IsAligned()) {
        *blockColumn = info.blockColumn;
        *blockWidth = info.block->m_width;
    }
}

Block * BlockMultipleAlignment::GetBlockBefore(const Block *block)
{
    Block *prevBlock = NULL;
    BlockList::iterator b, be = m_blocks.end();
    for (b=m_blocks.begin(); b!=be; ++b) {
        if (*b == block) break;
        prevBlock = b->GetPointer();
    }
    return prevBlock;
}

const Block * BlockMultipleAlignment::GetBlockBefore(const Block *block) const
{
    const Block *prevBlock = NULL;
    BlockList::const_iterator b, be = m_blocks.end();
    for (b=m_blocks.begin(); b!=be; ++b) {
        if (*b == block) break;
        prevBlock = b->GetPointer();
    }
    return prevBlock;
}

Block * BlockMultipleAlignment::GetBlockAfter(const Block *block)
{
    BlockList::iterator b, be = m_blocks.end();
    for (b=m_blocks.begin(); b!=be; ++b) {
        if (*b == block) {
            ++b;
            if (b == be) break;
            return *b;
        }
    }
    return NULL;
}

const Block * BlockMultipleAlignment::GetBlockAfter(const Block *block) const
{
    BlockList::const_iterator b, be = m_blocks.end();
    for (b=m_blocks.begin(); b!=be; ++b) {
        if (*b == block) {
            ++b;
            if (b == be) break;
            return *b;
        }
    }
    return NULL;
}

const UnalignedBlock * BlockMultipleAlignment::GetUnalignedBlockBefore(
    const UngappedAlignedBlock *aBlock) const
{
    const Block *prevBlock;
    if (aBlock)
        prevBlock = GetBlockBefore(aBlock);
    else
        prevBlock = m_blocks.back().GetPointer();
    return dynamic_cast<const UnalignedBlock*>(prevBlock);
}

void BlockMultipleAlignment::InsertBlockBefore(Block *newBlock, const Block *insertAt)
{
    BlockList::iterator b, be = m_blocks.end();
    for (b=m_blocks.begin(); b!=be; ++b) {
        if (*b == insertAt) {
            m_blocks.insert(b, CRef<Block>(newBlock));
            return;
        }
    }
    WARNING_MESSAGE("BlockMultipleAlignment::InsertBlockBefore() - couldn't find insertAt block");
}

void BlockMultipleAlignment::InsertBlockAfter(const Block *insertAt, Block *newBlock)
{
    BlockList::iterator b, be = m_blocks.end();
    for (b=m_blocks.begin(); b!=be; ++b) {
        if (*b == insertAt) {
            ++b;
            m_blocks.insert(b, CRef<Block>(newBlock));
            return;
        }
    }
    WARNING_MESSAGE("BlockMultipleAlignment::InsertBlockBefore() - couldn't find insertAt block");
}

void BlockMultipleAlignment::RemoveBlock(Block *block)
{
    BlockList::iterator b, be = m_blocks.end();
    for (b=m_blocks.begin(); b!=be; ++b) {
        if (*b == block) {
            m_blocks.erase(b);
            InitCache();
            return;
        }
    }
    WARNING_MESSAGE("BlockMultipleAlignment::RemoveBlock() - couldn't find block");
}

bool BlockMultipleAlignment::MoveBlockBoundary(unsigned int columnFrom, unsigned int columnTo)
{
    unsigned int blockColumn, blockWidth;
    GetAlignedBlockPosition(columnFrom, &blockColumn, &blockWidth);
    if (blockColumn == eUndefined || blockWidth == 0 || blockWidth == eUndefined) return false;

    TRACE_MESSAGE("trying to move block boundary from " << columnFrom << " to " << columnTo);

    const BlockInfo& info = m_blockMap[columnFrom];
    unsigned int row;
    int requestedShift = columnTo - columnFrom, actualShift = 0;
    const Block::Range *range;

    // shrink block from left
    if (blockColumn == 0 && requestedShift > 0 && requestedShift < (int) info.block->m_width) {
        actualShift = requestedShift;
        TRACE_MESSAGE("shrinking block from left");
        for (row=0; row<NRows(); ++row) {
            range = info.block->GetRangeOfRow(row);
            info.block->SetRangeOfRow(row, range->from + requestedShift, range->to);
        }
        info.block->m_width -= requestedShift;
        Block *prevBlock = GetBlockBefore(info.block);
        if (prevBlock && !prevBlock->IsAligned()) {
            for (row=0; row<NRows(); ++row) {
                range = prevBlock->GetRangeOfRow(row);
                prevBlock->SetRangeOfRow(row, range->from, range->to + requestedShift);
            }
            prevBlock->m_width += requestedShift;
        } else {
            Block *newUnalignedBlock = CreateNewUnalignedBlockBetween(prevBlock, info.block);
            if (newUnalignedBlock)
                InsertBlockBefore(newUnalignedBlock, info.block);
            TRACE_MESSAGE("added new unaligned block");
        }
    }

    // shrink block from right (requestedShift < 0)
    else if (blockColumn == info.block->m_width - 1 &&
                requestedShift < 0 && ((unsigned int) -requestedShift) < info.block->m_width) {
        actualShift = requestedShift;
        TRACE_MESSAGE("shrinking block from right");
        for (row=0; row<NRows(); ++row) {
            range = info.block->GetRangeOfRow(row);
            info.block->SetRangeOfRow(row, range->from, range->to + requestedShift);
        }
        info.block->m_width += requestedShift;
        Block *nextBlock = GetBlockAfter(info.block);
        if (nextBlock && !nextBlock->IsAligned()) {
            for (row=0; row<NRows(); ++row) {
                range = nextBlock->GetRangeOfRow(row);
                nextBlock->SetRangeOfRow(row, range->from + requestedShift, range->to);
            }
            nextBlock->m_width -= requestedShift;
        } else {
            Block *newUnalignedBlock = CreateNewUnalignedBlockBetween(info.block, nextBlock);
            if (newUnalignedBlock)
                InsertBlockAfter(info.block, newUnalignedBlock);
            TRACE_MESSAGE("added new unaligned block");
        }
    }

    // grow block to right
    else if (blockColumn == info.block->m_width - 1 && requestedShift > 0) {
        Block *nextBlock = GetBlockAfter(info.block);
        if (nextBlock && !nextBlock->IsAligned()) {
            int nRes;
            actualShift = requestedShift;
            for (row=0; row<NRows(); ++row) {
                range = nextBlock->GetRangeOfRow(row);
                nRes = range->to - range->from + 1;
                if (nRes < actualShift)
                    actualShift = nRes;
            }
            if (actualShift) {
                TRACE_MESSAGE("growing block to right");
                for (row=0; row<NRows(); ++row) {
                    range = info.block->GetRangeOfRow(row);
                    info.block->SetRangeOfRow(row, range->from, range->to + actualShift);
                    range = nextBlock->GetRangeOfRow(row);
                    nextBlock->SetRangeOfRow(row, range->from + actualShift, range->to);
                }
                info.block->m_width += actualShift;
                nextBlock->m_width -= actualShift;
                if (nextBlock->m_width == 0) {
                    RemoveBlock(nextBlock);
                    TRACE_MESSAGE("removed empty block");
                }
            }
        }
    }

    // grow block to left (requestedShift < 0)
    else if (blockColumn == 0 && requestedShift < 0) {
        Block *prevBlock = GetBlockBefore(info.block);
        if (prevBlock && !prevBlock->IsAligned()) {
            int nRes;
            actualShift = requestedShift;
            for (row=0; row<NRows(); ++row) {
                range = prevBlock->GetRangeOfRow(row);
                nRes = range->to - range->from + 1;
                if (nRes < -actualShift) actualShift = -nRes;
            }
            if (actualShift) {
                TRACE_MESSAGE("growing block to left");
                for (row=0; row<NRows(); ++row) {
                    range = info.block->GetRangeOfRow(row);
                    info.block->SetRangeOfRow(row, range->from + actualShift, range->to);
                    range = prevBlock->GetRangeOfRow(row);
                    prevBlock->SetRangeOfRow(row, range->from, range->to + actualShift);
                }
                info.block->m_width -= actualShift;
                prevBlock->m_width += actualShift;
                if (prevBlock->m_width == 0) {
                    RemoveBlock(prevBlock);
                    TRACE_MESSAGE("removed empty block");
                }
            }
        }
    }

    if (actualShift != 0) {
        UpdateBlockMap();
        return true;
    } else
        return false;
}

bool BlockMultipleAlignment::ShiftRow(unsigned int row, unsigned int fromAlignmentIndex, unsigned int toAlignmentIndex,
    eUnalignedJustification justification)
{
    if (fromAlignmentIndex == toAlignmentIndex)
        return false;

    Block
        *blockFrom = m_blockMap[fromAlignmentIndex].block,
        *blockTo = m_blockMap[toAlignmentIndex].block;

    // at least one end of the drag must be in an aligned block
    UngappedAlignedBlock *ABlock =
        dynamic_cast<UngappedAlignedBlock*>(blockFrom);
    if (ABlock) {
        if (blockTo != blockFrom && blockTo->IsAligned())
            return false;
    } else {
        ABlock = dynamic_cast<UngappedAlignedBlock*>(blockTo);
        if (!ABlock)
            return false;
    }

    // and the other must be in the same aligned block or an adjacent unaligned block
    UnalignedBlock
        *prevUABlock = dynamic_cast<UnalignedBlock*>(GetBlockBefore(ABlock)),
        *nextUABlock = dynamic_cast<UnalignedBlock*>(GetBlockAfter(ABlock));
    if (blockFrom != blockTo &&
        ((ABlock == blockFrom && prevUABlock != blockTo && nextUABlock != blockTo) ||
         (ABlock == blockTo && prevUABlock != blockFrom && nextUABlock != blockFrom)))