SequenceReader/SequenceReader.h

Go to the documentation of this file.

/*  Last edited: Apr 18 16:51 2002 (ac2) */

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

// 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  : SequenceReader
// File Name    : SequenceReader.h
// Language     : C++
// Module Author: Anthony J. Cox (ac2@sanger.ac.uk)

// Include guard:
#ifndef INCLUDED_SequenceReader
#define INCLUDED_SequenceReader

// Description:

// Includes:

// NB it is good practise for #include statements in header files to be
// replaced by forward declarations if at all possible
#include <GlobalDefinitions.h>
#include <string>
#include <iosfwd>
//class ostream;
class WordSequence;
class SequenceReader;
#include <HashTable.h>

// ### Class Declarations ###

class NumberOutOfRange : public SSAHAException
{
 public:
  NumberOutOfRange() :
    SSAHAException
    ("Requested sequence number exceeds number of sequences in file") {}
};

// Class Name : SequenceReaderMode
// Description: SequenceReaderMode and its subclasses encapsulate the 
// various policies for dealing with misread characters: ignore, replace with
// another character, report to monitoring stream. Subclasses define the
// actual behaviour, while SequenceReaderMode itself defines the common 
// interface 
class SequenceReaderMode
{
  public: 
  SequenceReaderMode( ostream& monStream = cout ):
  monitoringStream_( monStream )
  {}
  virtual ~SequenceReaderMode() {}
  // Check does three things:
  // Possibly: if thisChar is not valid, modify it to a valid base 
  // Possibly: print a message if thisChar is not valid
  // Always: return true if thisChar is to be processed, else false
  // Its exact behaviour depends on the subclass
  //  bool check( char& thisChar )
  //  {
  //    if (     (thisChar == 'A' ) || (thisChar == 'a') 
  //         || (thisChar == 'C' ) || (thisChar == 'c') 
  //          || (thisChar == 'G' ) || (thisChar == 'g') 
  //         || (thisChar == 'T' ) || (thisChar == 't') ) return true;
  //   return mismatch(thisChar);
  // }
  virtual bool mismatch( uchar& thisChar, Word& wordFlag ) const = 0;
  virtual SequenceReaderMode* clone( void ) = 0;
  protected:
  ostream& monitoringStream_;
}; // ~SequenceReaderMode

// Class Name : SequenceReaderModeIgnore
// Description: Simplest subclass of SequenceReaderMode - if a non-base 
// character is read, do nothing!
class SequenceReaderModeIgnore : public SequenceReaderMode
{
  public:
  SequenceReaderModeIgnore( ostream& monStream = cerr ):
  SequenceReaderMode( monStream )
  {}
  virtual bool mismatch( uchar& thisChar, Word& wordFlag ) const
  {
    DEBUG_L3("SequenceReaderModeIgnore::mismatch");
    return false;
  }
  virtual SequenceReaderMode* clone( void )
  {
    return new SequenceReaderModeIgnore(*this);
  } // ~clone

}; // ~SequenceReaderModeIgnore

// Class Name : SequenceReaderModeReport
// Description: Subclass of SequenceReaderMode - if a non-base 
// character is read, report to monitoring stream
class SequenceReaderModeReport : public SequenceReaderMode
{
  public:
  SequenceReaderModeReport( ostream& monStream = cout ):
  SequenceReaderMode( monStream )
  {}
  virtual bool mismatch( uchar& thisChar, Word& wordFlag ) const
  {
    DEBUG_L3("SequenceReaderModeReport::mismatch");
     
     monitoringStream_ << "Read unrecognized character (" 
                       << thisChar << ") from file" << endl;

     return false;
  } // ~mismatch
  virtual SequenceReaderMode* clone( void )
  {
    return new SequenceReaderModeReport(*this);
  } // ~clone


}; // ~SequenceReaderModeReport

// Class Name : SequenceReaderModeReplace
// Description: Subclass of SequenceReaderMode - if a non-base 
// character is read, silently replace it with a substitute
class SequenceReaderModeReplace : public SequenceReaderMode
{
  public:
  SequenceReaderModeReplace( uchar sub, ostream& monStream = cout ):
  SequenceReaderMode( monStream ), substitute_( sub )
  {}
  SequenceReaderModeReplace( const SequenceReaderModeReplace& rhs ) :
  substitute_( rhs.substitute_ ) {}

  virtual bool mismatch( uchar& thisChar, Word& wordFlag ) const
  {
    DEBUG_L3("SequenceReaderModeReplace::mismatch");
     if ( isgraph(thisChar) )
     {
       thisChar = substitute_;
       // Bug fix TC 14.9.00: now returns true, because thisChar is now valid
       // as the substitution has been done
       return true;
     }
     else return false;
  } // ~mismatch
  virtual SequenceReaderMode* clone( void )
  {
    return new SequenceReaderModeReplace(*this);
  } // ~clone
  protected:
  uchar substitute_;
}; // ~SequenceReaderModeReplace

// Class Name : SequenceReaderModeReportReplace
// Description: Subclass of SequenceReaderMode - if a non-base 
// character is read, replace it with a substitute and report to 
// monitoring stream
class SequenceReaderModeReportReplace : public SequenceReaderMode
{
  public:
  SequenceReaderModeReportReplace( uchar sub, ostream& monStream = cout ):
  SequenceReaderMode( monStream ), substitute_( sub )
  {}
  SequenceReaderModeReportReplace
  ( const SequenceReaderModeReportReplace& rhs ) :
  substitute_( rhs.substitute_ ) {}
  virtual bool mismatch( uchar& thisChar, Word& wordFlag ) const
  {
    DEBUG_L3("SequenceReaderModeReportReplace::mismatch");
     if ( isgraph(thisChar) )
     {
       monitoringStream_ << "Read unrecognized character (" 
                         << thisChar << ") from file, replacing with '" 
                         << substitute_ << "'.\n";
       thisChar = substitute_;
       // Bug fix TC 14.9.00: now returns true, because thisChar is now valid
       // as the substitution has been done
       return true;
     }
     else return false;
  } // ~mismatch
  virtual SequenceReaderMode* clone( void )
  {
    return new SequenceReaderModeReportReplace(*this);
  } // ~clone
  protected:
  uchar substitute_;
}; // ~SequenceReaderModeReportReplace

// Class Name : SequenceReaderModeFlagReplace
// Description: Subclass of SequenceReaderMode - if a non-base 
// character is read, silently replace it with a substitute and
// set wordFlag (whichSequenceEncoder will OR with the Word
// containing this letter 
class SequenceReaderModeFlagReplace : public SequenceReaderMode
{
  public:
  SequenceReaderModeFlagReplace( uchar sub, ostream& monStream = cout ):
  SequenceReaderMode( monStream ), substitute_( sub )
  {}
  SequenceReaderModeFlagReplace( const SequenceReaderModeFlagReplace& rhs ) :
  substitute_( rhs.substitute_ ) {}

  virtual bool mismatch( uchar& thisChar, Word& wordFlag ) const
  {
    DEBUG_L3("SequenceReaderModeFlagReplace::mismatch");
     if ( isgraph(thisChar) )
     {
       thisChar = substitute_;
       wordFlag = gCursedWord;
       return true;
     }
     else return false;
  } // ~mismatch
  virtual SequenceReaderMode* clone( void )
  {
    return new SequenceReaderModeFlagReplace(*this);
  } // ~clone
  protected:
  uchar substitute_;
}; // ~SequenceReaderModeFlagReplace





// -----------------------------

// Class Name : SequenceReaderPrinter
// Description: Purpose of this class and its subclasses is to provide an
// intuitive interface between sequence information and output streams.
// eg for SequenceReader myReader
// cout << myReader.getName(22);
// sends the name of sequence number 22 in myReader to standard output
// (assuming there are at least 22 sequences!)
// getName is a member object of SequenceReader that is an instance
// of SequenceReaderNamePrinter. By overloading the () operator we make
// it 'look like' a function.
class SequenceReaderPrinter
{
    public:
    // Constructor: called by the SequenceReader constructor 
    SequenceReaderPrinter( SequenceReader* inReader ):
      pReader_( inReader ) 
    {}

    SequenceReaderPrinter
    ( const SequenceReaderPrinter& rhs ) :
    pReader_( rhs.pReader_ ),
    seqNum_( rhs.seqNum_ ) {}


    SequenceReaderPrinter& operator()( SequenceNumber inSeqNum )
    { 
      seqNum_ = inSeqNum; return *this; 
    } // ~operator ()

    virtual void print( ostream& os ) = 0;

    friend ostream& operator<<(ostream& os, SequenceReaderPrinter& inPrinter )
    {
      inPrinter.print( os );
      return os;
    } // ~operator <<
  
    protected:
    SequenceReader* pReader_;
    SequenceNumber seqNum_;
}; // ~class SequenceReaderPrinter

// Class Name : SequenceReaderNamePrinter
// Description: Send the name of sequence seqNum_ to the ostream os
// by calling SequenceReader virtual member function printName
class SequenceReaderNamePrinter: public SequenceReaderPrinter
{
    public:
    SequenceReaderNamePrinter( SequenceReader* inReader ) :
    SequenceReaderPrinter( inReader ) {}
    
    virtual void print( ostream& os );

}; // ~class SequenceReaderNamePrinter

// Class Name : SequenceReaderSideInfoPrinter
// Description: Send the side info for sequence seqNum_ to the ostream os
// by calling SequenceReader virtual member function printSideInfo
class SequenceReaderSideInfoPrinter: public SequenceReaderPrinter
{
    public:
    SequenceReaderSideInfoPrinter( SequenceReader* inReader ) :
    SequenceReaderPrinter( inReader ) {}
    
    virtual void print( ostream& os );

}; // ~class SequenceReaderSideInfoPrinter

// Class Name : SequenceReaderSideInfoPrinter
// Description: Send the source file (ie the  raw ASCII) 
// for sequence seqNum_ to the ostream os by calling SequenceReader virtual 
// member function printSource
class SequenceReaderSourcePrinter: public SequenceReaderPrinter
{
    public:
    SequenceReaderSourcePrinter( SequenceReader* inReader ) :
    SequenceReaderPrinter( inReader ) {}
   
    virtual void print( ostream& os );

}; // ~SequenceReaderSourcePrinter


// ---------------------------
// Class Name : SourceReader
// Description: This is an abstract class that gives access to the original
// ASCII data that was used to produce the 2 bit per base sequence data
class SourceReader
{
 public:
  // if the position of the source buffer gets within 
  // resizeCacheThreshold_ of the current cache size then the
  // current cache size is multiplied by 1.5
  // The number of chars on a single line of fasta is thus 
  // effectively  limited to resizeCacheThreshold_ 
  enum
  { 
    sourceBufferSize_ = 20000, 
    resizeCacheThreshold_ = 5000, 
    nameBufferSize_ = 2000 
  };

  SourceReader( void ) : lastSourceSeqNum_(0) {}

  virtual ~SourceReader() {}
  // SeqInfo - data structure that stores the file number and position
  // of each indexed sequence 
  struct SeqIndexInfo
  {
    unsigned short fileNum;
    // NB this limits you to 2^32 ~= 4GB of seq in a single file
    // change seqPos to a std::streampos for more
    unsigned int seqPos;
  };

  // Function Name: extractSource
  // This extracts the source data for bases seqStart to seqEnd inclusive
  // of sequence seqNum and places it in source
  virtual void extractSource
    ( char** pSource, //vector<char>& source, 
    SequenceNumber seqNum,
    SequenceOffset seqStart,
    SequenceOffset seqEnd );

  // Function Name: extractSourceReverse
  // This extracts the source data using extractSource as above, then
  // reverse complements it.
  void extractSourceReverse
    ( char** pSource, //vector<char>& source, 
    SequenceNumber seqNum,
    SequenceOffset seqStart,
    SequenceOffset seqEnd );

  // Function Name: saveIndex
  // save the location of the start of each sequence in a file
  // sets up ofstream and calls the virtual saveIndexImp below
  void saveIndex( const string& fileName );

  // Function Name: extractToCache
  // Read in a sequence to lastSourceSeq_;
  void extractToCache( istream* pCurrentFile);

  // Function Name: saveIndexImp
  // Actually save the indexing data to disk. Implemented
  // for SequenceReaderFile and SequenceReaderMulti, not
  // for SourceReaderIndex
  virtual void saveIndexImp
  ( ostream& filesFile, 
    ostream& indexFile, 
    int& fileNumber );

 protected:
  int numCols_;
  // data to support SourceReader functionality
  vector<char> lastSourceSeq_;
  vector<char> reverseBuffer_;
  SequenceNumber lastSourceSeqNum_;
  char nameBuffer_[ nameBufferSize_ ];
};



// ---------------------------
// Class Name : SourceReaderIndex
// Description: Access source data using an index that gives entry points into
// a collection of files


class SourceReaderIndex : public SourceReader
{

  public:
  SourceReaderIndex( const string& fileName );
  virtual ~SourceReaderIndex();

  virtual void extractSource( char** pSource,//vector<char>& source, 
                     SequenceNumber seqNum,
                     SequenceOffset seqStart,
                     SequenceOffset seqEnd );

  const char* extractName( SequenceNumber seqNum );
  // NB no definition of saveIndex in this class
  SequenceNumber size( void ) const { return index_.size(); }
 protected:
  vector<string*> fileNames_;
  vector<SeqIndexInfo> index_;
  SequenceNumber currentFileNum_;
  char inputBuffer_[ sourceBufferSize_ ];
  string lastName_;
  // reading a new sequence changes lastNameSeqNum_ (since the name has to be
  // read to read the sequence_) but reading a new name does not change
  // lastSourceSeqNum_
  SequenceNumber lastNameSeqNum_;

  ifstream* pCurrentFile_;
  int numSeqs_;
};

class SequenceReaderState;

// -----------------------------

// Class Name : SequenceReader
// Description: This is an abstract class that specifies the interface via 
// via which the software will read the sequence data.
class SequenceReader : public SourceReader
{

  // PUBLIC MEMBER FUNCTIONS
  public:

  //  enum Constants{ notCalculatedYet_ = -9999 };

  // Constructors and Destructors

  // Function Name: Constructor
  // Arguments: ostream&
  SequenceReader( ostream& monitoringStream = cerr );

  // Function Name: Copy constructor
  // Arguments:
  // A copy constructor is required for all subclasses of SequenceReader.
  // This is because it is used by subclass SequenceReaderMulti (which forms
  // an aggregate of SequenceReader instances
  SequenceReader( const SequenceReader& rhs );


  // Function Name: Destructor
  // Arguments:
  virtual ~SequenceReader(); 
  // (NB destructor should be virtual if class is to be derived from)

  // Manipulator Functions
  virtual SequenceReader* clone( void ) = 0;


  // Function Name:
  // Arguments:
  // TYPE  NAME  IN/OUT COMMENT
  // Returns: TYPE COMMENT

  // Function Name: changeMode
  // Arguments: const SequenceReaderMode&
  // Makes a copy of mode and uses it to handle mismatch character reads
  virtual void changeMode( SequenceReaderMode* pMode ) = 0;


  // Accessor Functions
  // (NB all accessor functions should be 'const')

  // Function Name: rewind
  // Arguments: void
  // Returns:   void
  // Rewind to the start of the data file, so that getNextSequence will
  // return the first sequence in the file
  virtual void rewind( void ) = 0;

  // Function Name: findSequence
  // Arguments: SequenceNumber (in)
  // Returns:   void
  // Winds the input file stream to the start of sequence number seqNum. 
  // Returns false if seqNum exceeds the number of sequences in
  // the file.
  virtual bool findSequence( SequenceNumber seqNum ) = 0;

  // Function Name: getNextSequence
  // Arguments: WordSequence& (out), int (in)
  // Returns:   int
  // Read the next set of sequence information from the file and parse it
  // into WordSequence format. Returns -1 if there has been a problem with
  // reading the sequence, else returns the number of valid base pairs 
  // contained within the final word of the sequence.
  virtual int getNextSequence( WordSequence& nextSeq, int wordLength ) = 0;

  // Function Name: getNextSequence
  // Arguments: WordSequence& (out), HashTable& (in)
  // Returns:   int
  // Read the next set of sequence information from the file and parse it
  // into WordSequence format, using word length required by hashTable
  int getNextSequence( WordSequence& nextSeq, const HashTable& hashTable );

  // Function Name: getSequence
  // Arguments: WordSequence& (out), SequenceNumber (in), int (in)
  // Returns:   bool
  // Read the sequenceNumber-th set of sequence information from the file and 
  // parse it into WordSequence format
  virtual int getSequence
  ( WordSequence& nextSeq, SequenceNumber sequenceNumber, int wordLength ) = 0;

  // Function Name: getSequence
  // Arguments: WordSequence& (out), int (in), const HashTable& (in)
  // Returns:   bool
  // Read the sequenceNumber-th set of sequence information from the file and 
  // parse it into WordSequence format, getting word length from hashTable
  int getSequence
  ( WordSequence& nextSeq, 
    SequenceNumber sequenceNumber, 
    const HashTable& hashTable );

  // Function Name: getLastSequenceNumber
  // Arguments: void
  // Returns:   int
  // Returns the position in the data file of the last sequence read
  SequenceNumber getLastSequenceNumber( void ) const 
  { 
    return lastSequenceNumber_; 
  }

  // Function Name: areAllSequencesRead
  // Arguments: void
  // Returns:   bool
  // Returns true if the end of the file has been reached
  bool areAllSequencesRead( void ) const
  {
    return allSequencesRead_;
  }

  // Function Name: getNumSequencesInFile
  // Arguments: void
  // Returns:   SequenceNumber
  // Returns the number of sequences in the file (done by lazy 
  // initialization, i.e. will only be calculated if asked for).
  // NB current place in file will be lost.
  SequenceNumber getNumSequencesInFile( void )
  {
    if ( allSequencesRead_ == false )
    {
      numSequencesInFile_ = computeNumSequencesInFile();
      allSequencesRead_ = true;
    } 
    return numSequencesInFile_;
  } // ~getNumSequencesInFile( void )

  // Function Name: getLastSequenceName
  // Arguments: string& (out)
  // Returns:   void
  // Fills the string with the name of the last sequence read 
  virtual void getLastSequenceName( string& seqName ) const = 0;

  // Function Name: getBitsPerSymbol
  // Arguments: none
  // Returns:   int
  // Returns number of bits per symbol used in encoding
  virtual int getBitsPerSymbol ( void ) const = 0;

  // Function Name: getSourceDataType
  // Arguments: none
  // Returns:   SourceDataType
  // Returns type of data being encoded (protein or DNA)
  virtual SourceDataType getSourceDataType( void ) const = 0;

  SequenceReaderNamePrinter getName;
  SequenceReaderSideInfoPrinter getSideInfo;
  SequenceReaderSourcePrinter getSource;
 
  // Function Name: printName
  // Arguments: string& (out), SequenceNumber (in)
  // Returns:   void
  // Fills a string with the name of the requested sequence
  virtual bool printName( ostream& os, SequenceNumber seqNum ) = 0;

  // Function Name: printSideInfo
  // Arguments: string& (out), SequenceNumber (in)
  // Returns:   void
  // Fills a string with the name of the requested sequence
  virtual bool printSideInfo( ostream& os, SequenceNumber seqNum ) = 0;

  // Function Name: printSource
  // Arguments: string& (out), SequenceNumber (in)
  // Returns:   void
  // Fills a string with the name of the requested sequence
  virtual bool printSource( ostream& os, SequenceNumber seqNum ) = 0;

  // encodeBases task now done by SequenceEncoder::encode - TC 14.3.1
  // Function Name: encodeBases
  // Arguments: WordSequence& (out), const TYPE& (in), int (in), int(in)
  // Converts sequence data from character format into binary format and 
  // places it in seq. Making this function a template function means that
  // the same function can be used to read from character arrays or strings.
  //  template <class TYPE>
  //  void encodeBases
  // ( WordSequence& seq, const TYPE& data, int wordLength, int numChars, 
  //   int& basesInSequence = 0 );

  //void encodeBases
  //( WordSequence& seq, const string& data, int wordLength, int numChars,
  //  int& basesInSequence );
  //void encodeBases
  //( WordSequence& seq, const char* data, int wordLength, int numChars,
  // int& basesInSequence );

  // Functions to deal with state information

  // Function Name: saveState
  // Arguments: void
  // Returns:   SequenceReaderState*
  // saves the state (ie current file position) of a SequenceReader for future
  // restoration
  virtual SequenceReaderState* saveState( void ) const 
  { assert (1==0); return NULL;}

  // Function Name: restoreState
  // Arguments:   SequenceReaderState*
  // Returns:     void
  // restores the state (ie current file position) of a SequenceReader
  // then (NB!!) deletes *pState;
  virtual void restoreState( SequenceReaderState* pState ) 
  { assert (1==0);}




  // PROTECTED MEMBER FUNCTIONS 
  // (visible to this class and derived classes only)
  protected:

  // Function Name: computeNumSequencesInFile
  // Arguments: void
  // Returns:   int
  // Returns the number of sequences in the file - called by 
  // getNumSequencesInFile. NB this will lose the current place in the file.
  virtual SequenceNumber computeNumSequencesInFile( void ) = 0;


  // PRIVATE MEMBER FUNCTIONS
  // (visible to instances of this class only)
  
  private:
  SequenceReader& operator=(const SequenceReader&);   // NOT IMPLEMENTED

  // PROTECTED DATA:
  // (visible to instances of this class only)
  protected:

  // lastSequenceNumber_ is (surprise) the number of the last sequence that 
  // was read, and thus indicates the current position in the file. A value
  // of zero indicates we are at the beginning of the file
  SequenceNumber lastSequenceNumber_;     

  // allSequencesRead_ is initially set to false. Once the end of the file
  // has been reached for the first time it is set to true and 
  // numSequencesInFile_ is filled in.
  bool           allSequencesRead_;

  // numSequencesInFile_: we want to avoid calculating this unless we
  // absolutely have to, as it involves scanning all the way to the end of
  // the file, which can be slow for large files. 
  SequenceNumber numSequencesInFile_;
  ostream&  monitoringStream_; 
  //  SequenceReaderMode* pState_; now handled by SequenceEncoder

  // PRIVATE MEMBER DATA
  private:

}; // SequenceReader


// -----------------------------
// Class Name : SequenceReaderState
// Description: This preserves the state (ie position in a file) of
// a SequenceReader
class SequenceReaderState
{
 public:
  SequenceReaderState( SequenceNumber lsn ) :
    lastSequenceNumber_(lsn) {}
  virtual ~SequenceReaderState() {}
  const SequenceNumber lastSequenceNumber_; 
  // no point in making this private as it's const
};




// ### Function Declarations ###

// Name:
// Arguments:
// TYPE  NAME  IN/OUT COMMENT
// Returns: TYPE COMMENT

// End of include guard:
#endif

// End of file SequenceReader.h

---