cav_function.cpp

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

/*
 * ===========================================================================
 * PRODUCTION $Log: cav_function.cpp,v $
 * PRODUCTION Revision 1000.2  2004/06/01 19:41:19  gouriano
 * PRODUCTION PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.4
 * PRODUCTION
 * ===========================================================================
 */

/*  $Id: cav_function.cpp,v 1000.2 2004/06/01 19:41:19 gouriano Exp $
* ===========================================================================
*
*                            PUBLIC DOMAIN NOTICE
*               National Center for Biotechnology Information
*
*  This software/database is a "United States Government Work" under the
*  terms of the United States Copyright Act.  It was written as part of
*  the author's official duties as a United States Government employee and
*  thus cannot be copyrighted.  This software/database is freely available
*  to the public for use. The National Library of Medicine and the U.S.
*  Government have not placed any restriction on its use or reproduction.
*
*  Although all reasonable efforts have been taken to ensure the accuracy
*  and reliability of the software and data, the NLM and the U.S.
*  Government do not and cannot warrant the performance or results that
*  may be obtained by using this software or data. The NLM and the U.S.
*  Government disclaim all warranties, express or implied, including
*  warranties of performance, merchantability or fitness for any particular
*  purpose.
*
*  Please cite the author in any work or product based on this material.
*
* ===========================================================================
*
* Authors:  Paul Thiessen
*
* File Description:
*      C interfaced function body for cddalignview as function call
*
* ===========================================================================
*/

#include <ncbi_pch.hpp>
#include <corelib/ncbistl.hpp>
#include <corelib/ncbistre.hpp>
#include <corelib/ncbi_limits.h>
#include <corelib/ncbidiag.hpp>

#include <list>
#include <memory>

#include <objects/cdd/Cdd.hpp>
#include <objects/ncbimime/Ncbi_mime_asn1.hpp>
#include <objects/ncbimime/Biostruc_seqs.hpp>
#include <objects/ncbimime/Biostruc_align.hpp>
#include <objects/ncbimime/Biostruc_align_seq.hpp>

#include <objtools/cddalignview/cddalignview.h>
#include <objtools/cddalignview/cav_seqset.hpp>
#include <objtools/cddalignview/cav_alignset.hpp>
#include <objtools/cddalignview/cav_asnio.hpp>
#include <objtools/cddalignview/cav_alndisplay.hpp>


BEGIN_NCBI_SCOPE
USING_SCOPE(objects);

static EDiagSev defaultDiagPostLevel;

static int LoadASNFromIstream(CNcbiIstream& asnIstream,
    const SeqEntryList* *sequences, const SeqAnnotList* *alignments)
{
	*sequences = NULL;
	*alignments = NULL;

    // try to decide what ASN type this is, and if it's binary or ascii
    static const string
        asciiMimeFirstWord = "Ncbi-mime-asn1",
        asciiCDDFirstWord = "Cdd";
    bool isMime = false, isCDD = false, isBinary = true;

    string firstWord;
    asnIstream >> firstWord;
    if (firstWord == asciiMimeFirstWord) {
        isMime = true;
        isBinary = false;
    } else if (firstWord == asciiCDDFirstWord) {
        isCDD = true;
        isBinary = false;
    }

    // try to read the file as various ASN types (if it's not clear from the first ascii word).
    auto_ptr<SeqEntryList> newSequences(new SeqEntryList());
    auto_ptr<SeqAnnotList> newAlignments(new SeqAnnotList());
    bool readOK = false;
    string err;

    if (!isCDD) {
        ERR_POST(Info << "trying to read input as " <<
            ((isBinary) ? "binary" : "ascii") << " mime");
        CNcbi_mime_asn1 mime;
        SetDiagPostLevel(eDiag_Fatal); // ignore all but Fatal errors while reading data
        asnIstream.seekg(0);
        readOK = ReadASNFromIstream(asnIstream, mime, isBinary, err);
        SetDiagPostLevel(defaultDiagPostLevel);
        if (readOK) {
            // copy lists
            if (mime.IsStrucseqs()) {
                *newSequences = mime.GetStrucseqs().GetSequences();
                *newAlignments = mime.GetStrucseqs().GetSeqalign();
            } else if (mime.IsAlignstruc()) {
                *newSequences = mime.GetAlignstruc().GetSequences();
                *newAlignments = mime.GetAlignstruc().GetSeqalign();
            } else if (mime.IsAlignseq()) {
                *newSequences = mime.GetAlignseq().GetSequences();
                *newAlignments = mime.GetAlignseq().GetSeqalign();
            }
        } else {
            ERR_POST(Warning << "error: " << err);
        }
    }

    if (!readOK) {
        ERR_POST(Info << "trying to read input as " <<
            ((isBinary) ? "binary" : "ascii") << " cdd");
        CCdd cdd;
        SetDiagPostLevel(eDiag_Fatal); // ignore all but Fatal errors while reading data
        asnIstream.seekg(0);
        readOK = ReadASNFromIstream(asnIstream, cdd, isBinary, err);
        SetDiagPostLevel(defaultDiagPostLevel);
        if (readOK) {
            newSequences->resize(1);
            newSequences->front().Reset(&(cdd.SetSequences()));
            *newAlignments = cdd.GetSeqannot();   // copy the list
        } else {
            ERR_POST(Warning << "error: " << err);
        }
    }

    if (!readOK) {
        ERR_POST(Error << "Input is not a recognized data type (Ncbi-mime-asn1 or Cdd)");
        return CAV_ERROR_BAD_ASN;
    }

    *sequences = newSequences.release();
	*alignments = newAlignments.release();
    return CAV_SUCCESS;
}

// checks two things for each slave sequence: that all the residues of the sequence
// are present in the display, and that the aligned residues are in the right place
// wrt the master
static bool VerifyAlignmentData(const AlignmentSet *alignmentSet, const AlignmentDisplay *display)
{
    int alnLoc, masterLoc, slaveLoc, currentMasterLoc, currentSlaveLoc;
    char masterChar, slaveChar;
    const MasterSlaveAlignment *alignment;

    for (int i=0; i<alignmentSet->alignments.size(); ++i) {
        masterLoc = slaveLoc = -1;
        alignment = alignmentSet->alignments[i];

        for (alnLoc=0; alnLoc<display->GetWidth(); ++alnLoc) {

            // get and check characters
            masterChar = display->GetCharAt(alnLoc, 0);
            if (masterChar == '?') {
                ERR_POST(Error << "bad alignment coordinate: loc " << alnLoc << " row " << 0);
                return false;
            }
            slaveChar = display->GetCharAt(alnLoc, 1 + i);
            if (slaveChar == '?') {
                ERR_POST(Error << "bad alignment coordinate: loc " << alnLoc << " row " << (1+i));
                return false;
            }

            // advance seqLocs, check sequence string length and composition
            if (!IsGap(masterChar)) {
                ++masterLoc;
                if (i == 0) {   // only need to check master once
                    if (masterLoc >= alignment->master->sequenceString.size()) {
                        ERR_POST(Error << "master sequence too long at alnLoc " << alnLoc
                            << " row " << (i+1) << "masterLoc" << masterLoc);
                        return false;
                    } else if (toupper(masterChar) != toupper(alignment->master->sequenceString[masterLoc])) {
                        ERR_POST(Error << "master sequence mismatch at alnLoc " << alnLoc
                            << " row " << (i+1) << "masterLoc" << masterLoc);
                        return false;
                    }
                }
            }
            if (!IsGap(slaveChar)) {
                ++slaveLoc;
                if (slaveLoc >= alignment->slave->sequenceString.size()) {
                    ERR_POST(Error << "slave sequence too long at alnLoc " << alnLoc
                        << " row " << (i+1) << "slaveLoc" << slaveLoc);
                    return false;
                } else if (toupper(slaveChar) != toupper(alignment->slave->sequenceString[slaveLoc])) {
                    ERR_POST(Error << "slave sequence mismatch at alnLoc " << alnLoc
                        << " row " << (i+1) << "slaveLoc" << slaveLoc);
                    return false;
                }
            }
            currentMasterLoc = IsGap(masterChar) ? -1 : masterLoc;
            currentSlaveLoc = IsGap(slaveChar) ? -1 : slaveLoc;

            // check display characters, to see if they match alignment data
            if (IsGap(slaveChar) || IsUnaligned(slaveChar)) {
                if (currentMasterLoc >= 0 && alignment->masterToSlave[currentMasterLoc] != -1) {
                    ERR_POST(Error << "slave should be marked aligned at alnLoc " << alnLoc
                        << " row " << (i+1));
                    return false;
                }
            }
            if (IsAligned(slaveChar)) {
                if (!IsAligned(masterChar)) {
                    ERR_POST(Error <<" slave marked aligned but master unaligned at alnLoc " << alnLoc
                        << " row " << (i+1));
                    return false;
                }
                if (alignment->masterToSlave[currentMasterLoc] == -1) {
                    ERR_POST(Error << "slave incorrectly marked aligned at alnLoc " << alnLoc
                        << " row " << (i+1));
                    return false;
                }
                if (alignment->masterToSlave[currentMasterLoc] != currentSlaveLoc) {
                    ERR_POST(Error << "wrong slave residue aligned at alnLoc " << alnLoc
                        << " row " << (i+1));
                    return false;
                }
            }

            // converse: make sure alignment data is correctly reflected in display
            if (!IsGap(masterChar)) {
                if (alignment->masterToSlave[currentMasterLoc] == -1) {
                    if (IsAligned(slaveChar)) {