bl2seq.cpp

ncbi源码

/*
 * ===========================================================================
 * PRODUCTION $Log: bl2seq.cpp,v $
 * PRODUCTION Revision 1000.4  2004/06/01 18:05:35  gouriano
 * PRODUCTION PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.53
 * PRODUCTION
 * ===========================================================================
 */

/*  $Id: bl2seq.cpp,v 1000.4 2004/06/01 18:05:35 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.
 *
 * ===========================================================================
 *
 * Author:  Christiam Camacho
 *
 * ===========================================================================
 */

/// @file bl2seq.cpp
/// Implementation of CBl2Seq class.

#include <ncbi_pch.hpp>
#include <objmgr/util/sequence.hpp>
#include <objects/seqloc/Seq_loc.hpp>
#include <objects/seqalign/Seq_align_set.hpp>
#include <objects/seqalign/Seq_align.hpp>

#include <algo/blast/api/bl2seq.hpp>
#include <algo/blast/api/blast_options_handle.hpp>
#include <algo/blast/api/multiseq_src.hpp>
#include "blast_seqalign.hpp"
#include "blast_setup.hpp"

// NewBlast includes
#include <algo/blast/core/blast_def.h>
#include <algo/blast/core/blast_util.h>
#include <algo/blast/core/blast_setup.h>
#include <algo/blast/core/lookup_wrap.h>
#include <algo/blast/core/blast_engine.h>
#include <algo/blast/core/blast_traceback.h>

/** @addtogroup AlgoBlast
 *
 * @{
 */

BEGIN_NCBI_SCOPE
USING_SCOPE(objects);
BEGIN_SCOPE(blast)

CBl2Seq::CBl2Seq(const SSeqLoc& query, const SSeqLoc& subject, EProgram p)
    : mi_bQuerySetUpDone(false)
{
    TSeqLocVector queries;
    TSeqLocVector subjects;
    queries.push_back(query);
    subjects.push_back(subject);

    x_InitSeqs(queries, subjects);
    m_OptsHandle.Reset(CBlastOptionsFactory::Create(p));
}

CBl2Seq::CBl2Seq(const SSeqLoc& query, const SSeqLoc& subject,
                 CBlastOptionsHandle& opts)
    : mi_bQuerySetUpDone(false)
{
    TSeqLocVector queries;
    TSeqLocVector subjects;
    queries.push_back(query);
    subjects.push_back(subject);

    x_InitSeqs(queries, subjects);
    m_OptsHandle.Reset(&opts);
}

CBl2Seq::CBl2Seq(const SSeqLoc& query, const TSeqLocVector& subjects, 
                 EProgram p)
    : mi_bQuerySetUpDone(false)
{
    TSeqLocVector queries;
    queries.push_back(query);

    x_InitSeqs(queries, subjects);
    m_OptsHandle.Reset(CBlastOptionsFactory::Create(p));
}

CBl2Seq::CBl2Seq(const SSeqLoc& query, const TSeqLocVector& subjects, 
                 CBlastOptionsHandle& opts)
    : mi_bQuerySetUpDone(false)
{
    TSeqLocVector queries;
    queries.push_back(query);

    x_InitSeqs(queries, subjects);
    m_OptsHandle.Reset(&opts);
}

CBl2Seq::CBl2Seq(const TSeqLocVector& queries, const TSeqLocVector& subjects, 
                 EProgram p)
    : mi_bQuerySetUpDone(false)
{
    x_InitSeqs(queries, subjects);
    m_OptsHandle.Reset(CBlastOptionsFactory::Create(p));
}

CBl2Seq::CBl2Seq(const TSeqLocVector& queries, const TSeqLocVector& subjects, 
                 CBlastOptionsHandle& opts)
    : mi_bQuerySetUpDone(false)
{
    x_InitSeqs(queries, subjects);
    m_OptsHandle.Reset(&opts);
}

void CBl2Seq::x_InitSeqs(const TSeqLocVector& queries, 
                         const TSeqLocVector& subjs)
{
    m_tQueries = queries;
    m_tSubjects = subjs;
    mi_pScoreBlock = NULL;
    mi_pLookupTable = NULL;
    mi_pLookupSegments = NULL;
    mi_pResults = NULL;
    mi_pDiagnostics = NULL;
    mi_pSeqSrc = NULL;
}

CBl2Seq::~CBl2Seq()
{ 
    x_ResetQueryDs();
    x_ResetSubjectDs();
}

/// Resets query data structures
void
CBl2Seq::x_ResetQueryDs()
{
    mi_bQuerySetUpDone = false;
    // should be changed if derived classes are created
    mi_clsQueries.Reset(NULL);
    mi_clsQueryInfo.Reset(NULL);
    mi_pScoreBlock = BlastScoreBlkFree(mi_pScoreBlock);
    mi_pLookupTable = LookupTableWrapFree(mi_pLookupTable);
    mi_pLookupSegments = ListNodeFreeData(mi_pLookupSegments);
    // TODO: should clean filtered regions?
}

/// Resets subject data structures
void
CBl2Seq::x_ResetSubjectDs()
{
    // Clean up structures and results from any previous search
    mi_pSeqSrc = BlastSeqSrcFree(mi_pSeqSrc);
    mi_pResults = Blast_HSPResultsFree(mi_pResults);
    mi_pDiagnostics = Blast_DiagnosticsFree(mi_pDiagnostics);
    // TODO: Should clear class wrappers for internal parameters structures?
    //      -> destructors will be called for them
    //m_OptsHandle->SetDbSeqNum(0);  // FIXME: Really needed?
    //m_OptsHandle->SetDbLength(0);  // FIXME: Really needed?
}

TSeqAlignVector
CBl2Seq::Run() THROWS((CBlastException))
{
    //m_OptsHandle->GetOptions().DebugDumpText(cerr, "m_OptsHandle", 1);
    SetupSearch();
    m_OptsHandle->GetOptions().Validate();  // throws an exception on failure
    ScanDB();
    return x_Results2SeqAlign();
}


void 
CBl2Seq::SetupSearch()
{
    if ( !mi_bQuerySetUpDone ) {
        x_ResetQueryDs();
        SetupQueryInfo(m_tQueries, m_OptsHandle->GetOptions(), &mi_clsQueryInfo);
        SetupQueries(m_tQueries, m_OptsHandle->GetOptions(), mi_clsQueryInfo, 
                     &mi_clsQueries);

        // FIXME
        BlastMaskLoc* filter_mask = NULL;
        Blast_Message* blmsg = NULL;
        double scale_factor = 1.0;
        short st;

        st = BLAST_MainSetUp(m_OptsHandle->GetOptions().GetProgram(), 
                             m_OptsHandle->GetOptions().GetQueryOpts(),
                             m_OptsHandle->GetOptions().GetScoringOpts(),
                             m_OptsHandle->GetOptions().GetHitSaveOpts(),
                             mi_clsQueries, 
                             mi_clsQueryInfo, 
                             scale_factor,
                             &mi_pLookupSegments, 
                             &filter_mask, &mi_pScoreBlock, &blmsg);

        // Convert the BlastMaskLoc* into a CSeq_loc
        // TODO: Implement this! 
        //mi_vFilteredRegions = BLASTBlastMaskLoc2SeqLoc(filter_mask);
        BlastMaskLocFree(filter_mask); // FIXME, return seqlocs for formatter

        // TODO: Check that lookup_segments are not filtering the whole 
        // sequence (SSeqRange set to -1 -1)
        if (st != 0) {
            string msg = blmsg ? blmsg->message : "BLAST_MainSetUp failed";
            Blast_MessageFree(blmsg);
            NCBI_THROW(CBlastException, eInternal, msg.c_str());
        }
        Blast_MessageFree(blmsg);

        LookupTableWrapInit(mi_clsQueries, 
                            m_OptsHandle->GetOptions().GetLutOpts(),
                            mi_pLookupSegments, mi_pScoreBlock, 
                            &mi_pLookupTable, NULL);
        mi_bQuerySetUpDone = true;
    }

    x_ResetSubjectDs();

    mi_pSeqSrc = MultiSeqSrcInit(m_tSubjects, 
                                 m_OptsHandle->GetOptions().GetProgram());

    // Set the hitlist size to the total number of subject sequences, to 
    // make sure that no hits are discarded.
    m_OptsHandle->SetHitlistSize(m_tSubjects.size());
    m_OptsHandle->SetPrelimHitlistSize(m_tSubjects.size());
}

void 
CBl2Seq::ScanDB()