QueryManager/QueryManager.cpp

Go to the documentation of this file.

// #######################################################################

// SSAHA : Sequence Search and Alignment by Hashing Algorithm
// Version 3.2, released 1st March 2004
// Copyright (c) Genome Research 2002

// SSAHA is free software; you can redistribute it and/or modify 
// it under the terms of version 2 of the GNU General Public Licence
// as published by the Free Software Foundation.
 
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU General Public Licence for more details.
 
// You should have received a copy of the GNU General Public Licence
// along with this program; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
// or see the on-line version at http://www.gnu.org/copyleft/gpl.txt

// #######################################################################

// Module Name  : QueryManager
// File Name    : QueryManager.cpp
// Language     : C++
// Module Author: Anthony J. Cox (ac2@sanger.ac.uk)

// Description:

// Includes:

#include "QueryManager.h"
#include "SequenceReader.h"
#include "SequenceEncoder.h"
#include "HashTableGeneric.h"
#include "HashTableTranslated.h"
#include "MatchStoreUngapped.h"
#include "MatchStoreGapped.h"


// ### Function Definitions ###

// MatchAdder member function definitions

void MatchAdderImp::operator()(  SequenceNumber subjectNum,
                                     SequenceOffset numBases,
                                     SequenceOffset queryStart,
                                     SequenceOffset queryEnd,
                                     SequenceOffset subjectStart,
                                     SequenceOffset subjectEnd )
{
  //  if (subjectNum!=lastSubjectNum_)
  //   subjectTable_.getSequenceName(name_,subjectNum);
  //  lastSubjectNum_=subjectNum;
  pStore_->addMatch( // name_,
                     subjectNum, 
                     numBases, 
                     queryStart, 
                     queryEnd, 
                     subjectStart, 
                     subjectEnd, 
                     isQueryForward_, 
                     true );

} // ~void MatchAdderImp::operator()


MatchAdderCodonCodon::MatchAdderCodonCodon( HashTableTranslated& subjectTable ) :
MatchAdderImp( subjectTable ), subjectTable_( subjectTable ) {}

MatchAdderProteinCodon::MatchAdderProteinCodon( HashTableTranslated& subjectTable ) :
MatchAdderImp( subjectTable ), subjectTable_( subjectTable ) {}

MatchAdderCodonProtein::MatchAdderCodonProtein( HashTableGeneric& subjectTable ) :
MatchAdderImp( subjectTable ) {}


void MatchAdderCodonProtein::operator()(  SequenceNumber subjectNum,
                                     SequenceOffset numBases,
                                     SequenceOffset queryStart,
                                     SequenceOffset queryEnd,
                                     SequenceOffset subjectStart,
                                     SequenceOffset subjectEnd )
{
  //  if (subjectNum!=lastSubjectNum_)
  //    subjectTable_.getSequenceName(name_,subjectNum);
  //  lastSubjectNum_=subjectNum;
  pStore_->addMatch( // name_,
                     subjectNum, 
                     numBases, 
                     queryStart,
                     queryEnd, 
                     (subjectStart+2)/gNumReadingFrames, 
                     subjectEnd/gNumReadingFrames, 
                     isQueryForward_, true );
} // ~void MatchAdderCodonProtein::operator()

void MatchAdderCodonCodon::operator()(  SequenceNumber subjectNum,
                                     SequenceOffset numBases,
                                     SequenceOffset queryStart,
                                     SequenceOffset queryEnd,
                                     SequenceOffset subjectStart,
                                     SequenceOffset subjectEnd )
{
  if (subjectNum!=lastSubjectNum_)
  {
    // subjectTable_.getSequenceName(name_,subjectNum);
    size_ = subjectTable_.getSequenceSize(subjectNum);
    lastSubjectNum_=subjectNum;
  } // ~if
  pStore_->addMatch( // name_,
                     subjectNum,
                     numBases, 
                     queryStart,
                     queryEnd,
                     //              subjectStart, subjectEnd,
                     subjectTable_.isForward() ? subjectStart
                     : size_ - subjectEnd + 1,
                     subjectTable_.isForward() ? subjectEnd
                     : size_ - subjectStart + 1, 
                      isQueryForward_, subjectTable_.isForward() );
  

} // ~MatchAdderCodonCodon::operator()

void MatchAdderProteinCodon::operator()(  SequenceNumber subjectNum,
                                     SequenceOffset numBases,
                                     SequenceOffset queryStart,
                                     SequenceOffset queryEnd,
                                     SequenceOffset subjectStart,
                                     SequenceOffset subjectEnd )
{


  if (subjectNum!=lastSubjectNum_)
  {
    //    subjectTable_.getSequenceName(name_,subjectNum);
    size_ = subjectTable_.getSequenceSize(subjectNum);
    lastSubjectNum_=subjectNum;
  } // ~if

  //  cout << "MAPC: " << subjectNum << " " << queryStart << " " << queryEnd
  //   << " " << subjectStart << " " << subjectEnd << " - " << size_ << endl;

  pStore_->addMatch( // name_,
                     subjectNum, 
                     numBases/gNumReadingFrames, 
                     (queryStart+2)/gNumReadingFrames,  
                     queryEnd/gNumReadingFrames, 
                     //                              subjectStart,
                     //   subjectEnd,
                     subjectTable_.isForward() ? subjectStart
                     : size_ - subjectEnd + 1,
                     subjectTable_.isForward() ? subjectEnd
                     : size_ - subjectStart + 1, 
                     true,
                     subjectTable_.isForward() );
                     
} // ~MatchAdderProteinCodon::operator()



// MatchPolicy member function definitions

MatchPolicy::MatchPolicy( HashTableGeneric& subjectTable ) :
  subjectTable_( subjectTable ), 
  queryWordLength_( subjectTable.getWordLength() ) 
{}


MatchPolicyDNADNA::MatchPolicyDNADNA( HashTableGeneric& subjectTable ) :
MatchPolicy( subjectTable )
{
 addMatch_ = new MatchAdderImp(subjectTable_); 
}


void MatchPolicyDNADNA::operator()
  ( WordSequence& querySeqFwd, MatchStore& store, 
    MatchAlgorithm& findMatch )
{

  addMatch_->link(store);
  addMatch_->setQuerySize( ((querySeqFwd.size()-1) * queryWordLength_ )
      + querySeqFwd.getNumBasesInLast() );

  WordSequence querySeqRev;
  reverseComplement
    ( querySeqFwd, querySeqRev, 
      subjectTable_.getWordLength() ); 

  addMatch_->setQueryForward();

  findMatch( querySeqFwd, *addMatch_, subjectTable_ ); 

  addMatch_->setQueryReverse();

  findMatch( querySeqRev, *addMatch_, subjectTable_ );

} // ~QueryManager::doMatch


MatchPolicyProteinProtein::MatchPolicyProteinProtein
( HashTablePackedProtein& subjectTable ) :
subjectTable_( subjectTable ),
MatchPolicy( subjectTable )
{

  subjectTable_.setQueryProtein();
  addMatch_ = new MatchAdderImp(subjectTable_);
}

void MatchPolicyProteinProtein::operator()
  ( WordSequence& querySeqFwd, MatchStore& store, 
    MatchAlgorithm& findMatch )
{

  addMatch_->link(store);
  addMatch_->setQueryForward();

  findMatch( querySeqFwd, *addMatch_, subjectTable_ ); 

} // ~QueryManager::doMatch


MatchPolicyDNAProtein::MatchPolicyDNAProtein
( HashTablePackedProtein& subjectTable ) :
MatchPolicy( subjectTable ),
subjectTable_( subjectTable )
{
  queryWordLength_ = gMaxBasesPerWord;
  subjectTable_.setQueryTranslatedDNA();

  addMatch_ = new MatchAdderCodonProtein(subjectTable_);

} // ~MatchPolicyDNAProtein::MatchPolicyDNAProtein



void MatchPolicyDNAProtein::operator()
  ( WordSequence& querySeqFwd, MatchStore& store, 
    MatchAlgorithm& findMatch )
{

  addMatch_->link(store);
  addMatch_->setQueryForward();

  addMatch_->setQuerySize( ((querySeqFwd.size()-1) * gMaxBasesPerWord )
      + querySeqFwd.getNumBasesInLast() );

  WordSequence querySeqRev;
  reverseComplement
    ( querySeqFwd, querySeqRev, 
      queryWordLength_ ); 

  //  cout << "MPDP: doing fwd" << endl;

  addMatch_->setQueryForward();

  findMatch( querySeqFwd, *addMatch_, subjectTable_,
             gNumReadingFrames * subjectTable_.getWordLength() ); 

  //  cout << "MPDP: doing rev" << endl;

  addMatch_->setQueryReverse();

  findMatch( querySeqRev, *addMatch_, subjectTable_,
             gNumReadingFrames * subjectTable_.getWordLength() );


} // ~void MatchPolicyDNAProtein::operator()
 

// MatchPolicyProteinTranslated member functions

MatchPolicyProteinTranslated::MatchPolicyProteinTranslated
( HashTableTranslated& subjectTable ) :
MatchPolicy( subjectTable ),
subjectTable_( subjectTable )
{
  subjectTable_.setQueryProtein();
  addMatch_ = new MatchAdderProteinCodon(subjectTable_);
} // ~MatchPolicyProteinTranslated::MatchPolicyProteinTranslated


void MatchPolicyProteinTranslated::operator()
  ( WordSequence& querySeqFwd, MatchStore& store, 
    MatchAlgorithm& findMatch )
{

  addMatch_->link(store);
  addMatch_->setQueryForward();
  addMatch_->setQuerySize( ((querySeqFwd.size()-1) * gMaxBasesPerWord )
      + querySeqFwd.getNumBasesInLast() );

  // Necessary to copy querySeqFwd - otherwise it is shifted by the
  // first findMatch, which thus screws up the second findMatch
  WordSequence querySeqCopy(querySeqFwd);

  //  addMatch_->setSubjectForward();
  //  cout << "MPPT: fwd" << endl;
  subjectTable_.setForward();

  findMatch( querySeqCopy, *addMatch_, subjectTable_, 
             gNumReadingFrames * subjectTable_.getWordLength() );

  //  addMatch_->setSubjectReverse();
  // cout << "MPPT: rev" << endl;
  subjectTable_.setReverse();

  findMatch( querySeqFwd, *addMatch_, subjectTable_,
             gNumReadingFrames * subjectTable_.getWordLength() );


} // ~MatchPolicyProteinTranslated::operator()

// MatchPolicyDNATranslated member function defs

MatchPolicyDNATranslated::MatchPolicyDNATranslated
( HashTableTranslated& subjectTable ) :
MatchPolicy( subjectTable ),
subjectTable_( subjectTable )
{
   queryWordLength_ = gMaxBasesPerWord;
  subjectTable_.setQueryTranslatedDNA();
     addMatch_ = new MatchAdderCodonCodon(subjectTable_);
     //  addMatch_ = new MatchAdderImp(subjectTable_);
} // ~MatchPolicyDNATranslated::MatchPolicyDNATranslated




void MatchPolicyDNATranslated::operator()
  ( WordSequence& querySeqFwd, MatchStore& store, 
    MatchAlgorithm& findMatch )
{

  addMatch_->link(store);

  addMatch_->setQuerySize( ((querySeqFwd.size()-1) * gMaxBasesPerWord )
      + querySeqFwd.getNumBasesInLast() );

  WordSequence revSeq, translatedQuery;

  reverseComplement( querySeqFwd, revSeq, gMaxBasesPerWord );

  addMatch_->setQueryForward();

  //  addMatch_->setSubjectForward();
  subjectTable_.setForward();

  findMatch( querySeqFwd, *addMatch_, subjectTable_, 
             gNumReadingFrames * subjectTable_.getWordLength() );

  //  addMatch_->setSubjectReverse();
  subjectTable_.setReverse();

  findMatch( querySeqFwd, *addMatch_, subjectTable_, 
             gNumReadingFrames * subjectTable_.getWordLength() );

  addMatch_->setQueryReverse();

  //  addMatch_->setSubjectForward();
  subjectTable_.setForward();

  findMatch( revSeq, *addMatch_, subjectTable_, 
             gNumReadingFrames * subjectTable_.getWordLength() );

  //  addMatch_->setSubjectReverse();
  subjectTable_.setReverse();

  findMatch( revSeq, *addMatch_, subjectTable_, 
             gNumReadingFrames * subjectTable_.getWordLength() );



} // ~MatchPolicyDNATranslated::operator()

// QueryManager function definitions

// Name:
// Arguments:
// TYPE  NAME  IN/OUT COMMENT
// Returns: TYPE COMMENT
  QueryManager::QueryManager
  ( SequenceReader& querySeqs,  
    HashTableGeneric& subjectSeqs, ostream& monitoringStream ) :
   queryReader_( querySeqs ), 
   subjectTable_(  subjectSeqs ), 
   monitoringStream_( monitoringStream )
   {
     monitoringStream_ << "constructing QueryManager\n";

     HashTableTranslated* pTransTable
       (dynamic_cast<HashTableTranslated*>(&subjectSeqs));

     if ( pTransTable != NULL )
     {
       if (queryReader_.getBitsPerSymbol()==gResidueBits )
       {
         monitoringStream_ 
         << "Info: running protein query against translated DNA hash table.\n";
         policy_ = new MatchPolicyProteinTranslated(*pTransTable);
       }         
       else if (queryReader_.getBitsPerSymbol()==gBaseBits )
       {
         monitoringStream_ 
         << "Info: running DNA query against translated DNA hash table.\n";
         policy_ = new MatchPolicyDNATranslated(*pTransTable);
       }
       else assert (1==0);
     }
     else // not a translated table
     {
       HashTablePackedProtein* pProteinTable
         (dynamic_cast<HashTablePackedProtein*>(&subjectSeqs));
       if (pProteinTable!=NULL)
       {

         if (queryReader_.getBitsPerSymbol()==gResidueBits)
         {
           monitoringStream_ 
             << "Info: running protein query against protein hash table."
             << endl;
           policy_ = new MatchPolicyProteinProtein(*pProteinTable);
         } // ~if
         else if (queryReader_.getBitsPerSymbol()==gBaseBits)
         {
           monitoringStream_ 
             << "Info: running DNA query against protein hash table."
             << endl;
           policy_ = new MatchPolicyDNAProtein(*pProteinTable);
         } // ~else if
         else assert (1==0);  
       } // ~if  
       else if (subjectTable_.getBitsPerSymbol()==gBaseBits )
       {
         if (queryReader_.getBitsPerSymbol()==gBaseBits)
         {
           monitoringStream_ 
           << "Info: running DNA query against DNA hash table."
           << endl;
           policy_ = new MatchPolicyDNADNA(subjectSeqs);
         } // ~if
         else
         {
           monitoringStream_ 
           << "Error: can't run protein query against DNA hash table!\n";
           throw SSAHAException
             ("Can't run protein query against DNA hash table!");
         } // ~else
         
       } // ~else if

     } // ~else

   } // ~constructor

   QueryManager::~QueryManager()
   {
     monitoringStream_ << "destructing QueryManager\n";
     delete policy_;
   } // ~constructor

  void QueryManager::doQuery
  //  ( MatchStore& store, 
  ( MatchAlgorithm& match, 
    MatchTask& task,
    int queryStart, 
    int queryEnd )
  {
    int numBasesInLast(0);
    int wordLength( (*policy_).getWordLength() );

    WordSequence querySeqFwd;

    queryReader_.rewind();

    if ( queryStart == 1 ) 
    {
      numBasesInLast = queryReader_.getNextSequence
                         ( querySeqFwd, wordLength );
    } // ~if 
    else if ( queryStart > 1 ) 
    {
      numBasesInLast = queryReader_.getSequence
      ( querySeqFwd, queryStart, wordLength );
    } // ~else if
    else throw SSAHAException("Invalid query start value!");

    if ( numBasesInLast == -1 )
    {
      monitoringStream_ << "Info: requested sequence start (" << queryStart
      << ") exceeds number of sequences in query database.\n";
      return;
    } // ~if

    do
    {

      string queryName;
      queryReader_.getLastSequenceName( queryName );
      int queryNum = queryReader_.getLastSequenceNumber();
      int queryBases = ( querySeqFwd.size() - 1 )*wordLength 
        + querySeqFwd.getNumBasesInLast();

      MatchStoreImp store
      ( queryName, 
        queryNum, 
        queryBases,  
        subjectTable_.getNameReader() );

      //      cout << "QM: doing " << queryNum << endl;
      (*policy_)( querySeqFwd, store, match );

      store.setup();

      task( store );

      if ( queryReader_.getLastSequenceNumber() == queryEnd ) break;

      // clear the query sequence ready to read in next query
      querySeqFwd.clear();


      // read in next query
      numBasesInLast = queryReader_.getNextSequence
                         ( querySeqFwd, wordLength );
    } // ~while
    while ( numBasesInLast != -1 );

    if (    ( queryReader_.getLastSequenceNumber() < queryEnd  )
         && ( queryEnd != - 1 ) )
    {
      monitoringStream_ << "Info: requested final sequence (" << queryEnd
      << ") exceeded number of\nsequences in query database ("
      << queryReader_.getLastSequenceNumber() << ").\n";
    } // ~if
    
    return;

  } // ~QueryManager::doQuery




// End of file QueryManager.cpp

---