alignment_manager.cpp

ncbi源码

            return NULL;
        }
    }
    BlockMultipleAlignment *multipleAlignment = new BlockMultipleAlignment(sequenceList, this);

    // each block is a continuous region on the master, over which each master
    // residue is aligned to a residue of each slave, and where there are no
    // insertions relative to the master in any of the slaves
    int masterFrom = 0, masterTo, row;
    UngappedAlignedBlock *newBlock;

    while (masterFrom < multipleAlignment->GetMaster()->Length()) {

        // look for first all-aligned residue
        if (!AlignedToAllSlaves(masterFrom, alignments)) {
            ++masterFrom;
            continue;
        }

        // find all next continuous all-aligned residues, but checking for
        // block boundaries from the original master-slave pairs, so that
        // blocks don't get merged
        for (masterTo=masterFrom+1;
                masterTo < multipleAlignment->GetMaster()->Length() &&
                AlignedToAllSlaves(masterTo, alignments) &&
                NoSlaveInsertionsBetween(masterFrom, masterTo, alignments) &&
                NoBlockBoundariesBetween(masterFrom, masterTo, alignments);
             ++masterTo) ;
        --masterTo; // after loop, masterTo = first residue past block

        // create new block with ranges from master and all slaves
        newBlock = new UngappedAlignedBlock(multipleAlignment);
        newBlock->SetRangeOfRow(0, masterFrom, masterTo);
        newBlock->width = masterTo - masterFrom + 1;

        //TESTMSG("masterFrom " << masterFrom+1 << ", masterTo " << masterTo+1);
        for (a=alignments.begin(), row=1; a!=ae; ++a, ++row) {
            newBlock->SetRangeOfRow(row,
                (*a)->masterToSlave[masterFrom],
                (*a)->masterToSlave[masterTo]);
            //TESTMSG("slave->from " << b->from+1 << ", slave->to " << b->to+1);
        }

        // copy new block into alignment
        multipleAlignment->AddAlignedBlockAtEnd(newBlock);

        // start looking for next block
        masterFrom = masterTo + 1;
    }

    if (!multipleAlignment->AddUnalignedBlocks() ||
        !multipleAlignment->UpdateBlockMapAndColors()) {
        ERRORMSG("AlignmentManager::CreateMultipleFromPairwiseWithIBM() - "
            "error finalizing alignment");
        return NULL;
    }

    return multipleAlignment;
}

static void GetAlignedResidueIndexes(
    BlockMultipleAlignment::UngappedAlignedBlockList::const_iterator& b,
    BlockMultipleAlignment::UngappedAlignedBlockList::const_iterator& be,
    int row, int *seqIndexes)
{
    int i = 0, c;
    const Block::Range *range;
    for (; b!=be; ++b) {
        range = (*b)->GetRangeOfRow(row);
        for (c=0; c<(*b)->width; ++c) {
            seqIndexes[i++] = range->from + c;
        }
    }
}

void AlignmentManager::RealignAllSlaveStructures(void) const
{
    const BlockMultipleAlignment *multiple = GetCurrentMultipleAlignment();
    if (!multiple) return;
    BlockMultipleAlignment::UngappedAlignedBlockList blocks;
    multiple->GetUngappedAlignedBlocks(&blocks);
    if (blocks.size() == 0) {
        WARNINGMSG("Can't realign slaves with no aligned residues!");
        return;
    }
    BlockMultipleAlignment::UngappedAlignedBlockList::const_iterator b, be = blocks.end();
    int nResidues = 0;
    for (b=blocks.begin(); b!=be; ++b) nResidues += (*b)->width;
    if (nResidues <= 2) {
        WARNINGMSG("Can't realign slaves with < 3 aligned residues!");
        return;
    }

    const Sequence *masterSeq = multiple->GetSequenceOfRow(0), *slaveSeq;
    const Molecule *masterMol, *slaveMol;
    if (!masterSeq || !(masterMol = masterSeq->molecule)) {
        WARNINGMSG("Can't realign slaves to non-structured master!");
        return;
    }

    int *masterSeqIndexes = new int[nResidues], *slaveSeqIndexes = new int[nResidues];
    b = blocks.begin();
    GetAlignedResidueIndexes(b, be, 0, masterSeqIndexes);

    double *weights = new double[nResidues];

    const StructureObject *slaveObj;

    typedef const Vector * CVP;
    CVP *masterCoords = new CVP[nResidues], *slaveCoords = new CVP[nResidues];
    if (!masterMol->GetAlphaCoords(nResidues, masterSeqIndexes, masterCoords)) {
        WARNINGMSG("Can't get master alpha coords");
    } else if (masterSeq->GetOrSetMMDBLink() == MoleculeIdentifier::VALUE_NOT_SET) {
        WARNINGMSG("Don't know master MMDB ID");
    } else {

        masterMol->parentSet->InitStructureAlignments(masterSeq->identifier->mmdbID);

        int nStructureAlignments = 0;
        for (int i=1; i<multiple->NRows(); ++i) {
            slaveSeq = multiple->GetSequenceOfRow(i);
            if (!slaveSeq || !(slaveMol = slaveSeq->molecule)) continue;

            b = blocks.begin();
            GetAlignedResidueIndexes(b, be, i, slaveSeqIndexes);
            if (slaveMol->GetAlphaCoords(nResidues, slaveSeqIndexes, slaveCoords) < 3) {
                ERRORMSG("can't realign slave " << slaveSeq->identifier->pdbID << ", not enough coordinates in aligned region");
                continue;
            }

            if (!slaveMol->GetParentOfType(&slaveObj)) continue;

            // if any Vector* is NULL, make sure that weight is 0 so the pointer won't be accessed
            for (int j=0; j<nResidues; ++j) {
                if (!masterCoords[j] || !slaveCoords[j])
                    weights[j] = 0.0;
                else
                    weights[j] = 1.0; // for now, just use flat weighting
            }

            INFOMSG("realigning slave " << slaveSeq->identifier->pdbID << " against master " << masterSeq->identifier->pdbID);
            (const_cast<StructureObject*>(slaveObj))->
                RealignStructure(nResidues, masterCoords, slaveCoords, weights, i);
            ++nStructureAlignments;
        }

        // if no structure alignments, remove the list entirely
        if (nStructureAlignments == 0) masterMol->parentSet->RemoveStructureAlignments();
    }

    delete[] masterSeqIndexes;
    delete[] slaveSeqIndexes;
    delete[] masterCoords;
    delete[] slaveCoords;
    delete[] weights;
    return;
}

void AlignmentManager::GetAlignmentSetSlaveSequences(vector < const Sequence * > *sequences) const
{
    sequences->resize(alignmentSet->alignments.size());

    AlignmentSet::AlignmentList::const_iterator a, ae = alignmentSet->alignments.end();
    int i = 0;
    for (a=alignmentSet->alignments.begin(); a!=ae; ++a, ++i) {
        (*sequences)[i] = (*a)->slave;
    }
}

void AlignmentManager::GetAlignmentSetSlaveVisibilities(vector < bool > *visibilities) const
{
    if (slavesVisible.size() != alignmentSet->alignments.size()) // can happen if row is added/deleted
        slavesVisible.resize(alignmentSet->alignments.size(), true);

    // copy visibility list
    *visibilities = slavesVisible;
}

void AlignmentManager::ShowHideCallbackFunction(const vector < bool >& itemsEnabled)
{
    if (itemsEnabled.size() != slavesVisible.size() ||
        itemsEnabled.size() != alignmentSet->alignments.size()) {
        ERRORMSG("AlignmentManager::ShowHideCallbackFunction() - wrong size list");
        return;
    }

    slavesVisible = itemsEnabled;
    NewMultipleWithRows(slavesVisible);

    AlignmentSet::AlignmentList::const_iterator
        a = alignmentSet->alignments.begin(), ae = alignmentSet->alignments.end();
    const StructureObject *object;

    if ((*a)->master->molecule) {
        // Show() redraws whole StructureObject only if necessary
        if ((*a)->master->molecule->GetParentOfType(&object))
            object->parentSet->showHideManager->Show(object, true);
        // always redraw aligned molecule, in case alignment colors change
        GlobalMessenger()->PostRedrawMolecule((*a)->master->molecule);
    }
    for (int i=0; a!=ae; ++a, ++i) {
        if ((*a)->slave->molecule) {
            if ((*a)->slave->molecule->GetParentOfType(&object))
                object->parentSet->showHideManager->Show(object, slavesVisible[i]);
            GlobalMessenger()->PostRedrawMolecule((*a)->slave->molecule);
        }
    }

    // do necessary redraws + show/hides: sequences + chains in the alignment
    sequenceViewer->Refresh();
    GlobalMessenger()->PostRedrawAllSequenceViewers();
    GlobalMessenger()->UnPostRedrawSequenceViewer(sequenceViewer);  // Refresh() does this already
}

void AlignmentManager::NewMultipleWithRows(const vector < bool >& visibilities)
{
    if (visibilities.size() != alignmentSet->alignments.size()) {
        ERRORMSG("AlignmentManager::NewMultipleWithRows() - wrong size visibility vector");
        return;
    }

    // make a multiple from all visible rows
    PairwiseAlignmentList alignments;
    AlignmentSet::AlignmentList::const_iterator a, ae=alignmentSet->alignments.end();
    int i = 0;
    for (a=alignmentSet->alignments.begin(); a!=ae; ++a, ++i) {
        if (visibilities[i])
            alignments.push_back(*a);
    }

    // sequenceViewer will own the resulting alignment
    sequenceViewer->DisplayAlignment(CreateMultipleFromPairwiseWithIBM(alignments));
}

bool AlignmentManager::IsAligned(const Sequence *sequence, int seqIndex) const
{
    const BlockMultipleAlignment *currentAlignment = GetCurrentMultipleAlignment();
    if (currentAlignment)
        return currentAlignment->IsAligned(sequence, seqIndex);
    else
        return false;
}

bool AlignmentManager::IsInAlignment(const Sequence *sequence) const
{
    if (!sequence) return false;
    const BlockMultipleAlignment *currentAlignment = GetCurrentMultipleAlignment();
    if (currentAlignment) {
        for (int i=0; i<currentAlignment->sequences->size(); ++i) {
            if ((*(currentAlignment->sequences))[i] == sequence)
                return true;
        }
    }
    return false;
}

const Vector * AlignmentManager::GetAlignmentColor(const Sequence *sequence, int seqIndex,
    StyleSettings::eColorScheme colorScheme) const
{
    const BlockMultipleAlignment *currentAlignment = GetCurrentMultipleAlignment();
    if (currentAlignment)
        return currentAlignment->GetAlignmentColor(sequence, seqIndex, colorScheme);
    else
        return NULL;
}

void AlignmentManager::ShowSequenceViewer(void) const
{
    sequenceViewer->CreateSequenceWindow(true);
}

void AlignmentManager::ShowUpdateWindow(void) const
{
    updateViewer->CreateUpdateWindow();
}

void AlignmentManager::RealignSlaveSequences(
    BlockMultipleAlignment *multiple, const vector < int >& slavesToRealign)
{
    if (!multiple || sequenceViewer->GetCurrentAlignments().size() == 0 ||
            multiple != sequenceViewer->GetCurrentAlignments().front()) {
        ERRORMSG("AlignmentManager::RealignSlaveSequences() - wrong multiple alignment");
        return;
    }
    if (slavesToRealign.size() == 0) return;

    // create alignments for each master/slave pair, then update displays
    UpdateViewer::AlignmentList alignments;
    TRACEMSG("extracting rows");
    if (multiple->ExtractRows(slavesToRealign, &alignments)) {
        TRACEMSG("recreating display");
        sequenceViewer->GetCurrentDisplay()->RowsRemoved(slavesToRealign, multiple);
        TRACEMSG("adding to update window");
        SetDiagPostLevel(eDiag_Warning);    // otherwise, info messages take a long time if lots of rows
        updateViewer->AddAlignments(alignments);
        SetDiagPostLevel(eDiag_Info);
        TRACEMSG("done");
        updateViewer->CreateUpdateWindow();
    }
}

void AlignmentManager::ThreadUpdate(const ThreaderOptions& options, BlockMultipleAlignment *single)
{
    const ViewerBase::AlignmentList& currentAlignments = sequenceViewer->GetCurrentAlignments();
    if (currentAlignments.size() == 0) return;

    // make sure the editor is on in the sequenceViewer if merge is selected
    if (!sequenceViewer->EditorIsOn() && options.mergeAfterEachSequence) {
        ERRORMSG("Can only merge updates when editing is enabled in the sequence window");
        return;
    }

    // run the threader on the given alignment
    UpdateViewer::AlignmentList singleList, replacedList;