C++ Class: SamFile

Reading/Writing SAM/BAM Files In Your Program

The SamFile class allows a user to easily read/write a SAM/BAM file.

The SamFile class contains additional functionality that allows a user to read specific sections of sorted & indexed BAM files. In order take advantage of this capability, the index file must be read prior to setting the read section. This logic saves the time of having to read the entire file and takes advantage of the seeking capability of BGZF files.

Future Enhancements: Add the ability to read alignments that match a given start, end position for a specific reference sequence.

This class is part of libbam.

Class Methods

Method Name	Description
bool SamFile::IsEOF()	bool: true if the end of file has been reached, false if not. Be careful using this method when you are only reading a specific section since you may reach the end of your section without hitting the end of the file
bool SamFile::OpenForRead(const char* filename)	Opens the specified file for reading. Determines if it is a BAM/SAM file by reading the beginning of the file. Returns true if successfully opened reading, false if not.
bool SamFile::OpenForWrite(const char * filename)	bool: true if successfully opened, false if not. Opens as BAM file if the specified filename ends in .bam. Otherwise it is opened as a SAM file. Returns true if successfully opened for writing, false if not.
bool SamFile::ReadBamIndex(const char * filename)	bool: true if the bam index file was successfully read, false if not. Reads the specified bam index file. It must be read prior to setting a read section, for seeking and reading portions of a bam file.
void SamFile::Close()	Close the file if there is one open.
bool SamFile::ReadHeader(SamFileHeader& header)	Reads the header section from the file and stores it in the passed in header. Returns true if successfully read, false if not.
bool SamFile::WriteHeader(SamFileHeader& header)	Writes the specified header into the file. Returns true if successfully written, false if not.
bool SamFile::ReadRecord(SamFileHeader& header, SamRecord& record)	Reads the next record from the file and stores it in the passed in record. If it is an indexed BAM file and SetReadSection was called, only alignments in the section specified by SetReadSection are read. If they all have already been read, this method returns false. Validates that the record is sorted according to the value set by setSortedValidation. No sorting validation is done if specified to be unsorted, or setSortedValidation was never called. Returns false if the record was not successfully read or not properly sorted. Returns true if successfully read and properly sorted.
bool SamFile::WriteRecord(SamFileHeader& header, SamRecord& record)	Writes the specified record into the file. Validates that the record is sorted according to the value set by setSortedValidation. No sorting validation is done if specified to be unsorted, or setSortedValidation was never called. Returns false and does not write the record if the record was not properly sorted. Returns false if the record was not properly sorted or not successfully written. Returns true if properly sorted and successfully written.
void SamFile::setSortedValidation(SortedType sortType)<\code>	Set the flag to validate that the file is sorted as it is read/written. Must be called after the file has been opened. sortType specifies the type of sort to be checked for.
uint32_t SamFile::GetCurrentRecordCount()	Return the number of records that have been read/written so far.
SamStatus::Status SamFile::GetStatus()	Get the status result of the last status reporting method call.
const char* SamFile::GetStatusMessage()	Get the Status Message of the last call that sets status.
DEPRECATED: SamStatus::Status SamFile::GetFailure()	Get the type of failure that occurred on a method failure.
bool SamFile::SetReadSection(int32_t refID)	Tell the class which reference ID should be read from the BAM file. This is the index into the BAM Index list of reference information: 0 - #references. The records for that reference id will be retrieved on each ReadRecord call. When all records have been retrieved for the specified reference id, ReadRecord will return false until a new read section is set. Pass in -1 in order to read the section of the bam file not associated with any reference ID. Returns true if the read section was successfully set, false if not. False is returned if the BAM Index File has not yet been read or if a BAM file is not open for reading.

Class Enums

enum SortedType
Enum Value	Description
UNSORTED	Do not do any sorting validation.
FLAG	Validate that the file is sorted by the type specified in the SO Tag of the file's header.
COORDINATE	Validate that the file is sorted by Coordinate.
QUERY_NAME	Validate that the file is sorted by Query Name.

Usage Examples

Simple Read/Write SAM/BAM files

The following example reads in a sam/bam file and writes it out as a sam/bam file. The file format of the input sam/bam is determined by the SamFile class based on reading the type from the file. The file format of the output sam/bam file is determined by the SamFile class based on the extension of the output file. A ".bam" extension indicates a BAM file. All other extensions indicate SAM files.

#include "SamFile.h"
#include "SamValidation.h"

int main(int argc, char ** argv)
{
   // Check for the appropriate number of arguments.
   if(argc != 3)
   {
      printf("./bam <inputFile> <outputFile.sam/bam>\n");
      exit(-1);
   }

   // Open the input file for reading.
   SamFile samIn;
   if(!samIn.OpenForRead(argv[1]))
   {
      fprintf(stderr, "%s\n", samIn.GetStatusMessage());
      return(samIn.GetStatus());
   }

   // Open the output file for writing.
   SamFile samOut;
   if(!samOut.OpenForWrite(argv[2]))
   {
      fprintf(stderr, "%s\n", samOut.GetStatusMessage());
      return(samOut.GetStatus());
   }

   // Read the sam header.
   SamFileHeader samHeader;
   if(!samIn.ReadHeader(samHeader))
   {
      fprintf(stderr, "%s\n", samIn.GetStatusMessage());
      return(samIn.GetStatus());
   }

   // Write the sam header.
   if(!samOut.WriteHeader(samHeader))
   {
      fprintf(stderr, "%s\n", samOut.GetStatusMessage());
      return(samOut.GetStatus());     
   }

   SamRecord samRecord;
   SamValidationErrors samInvalidErrors;

   // Set writeStatus to success.  It will be changed
   // to the failure reason if any of the writes fail.
   SamStatus::Status writeStatus = SamStatus::SUCCESS;

   // Keep reading records until ReadRecord returns false.
   while(samIn.ReadRecord(samHeader, samRecord))
   {
      // Successfully read a record from the file, so write it.
      if(!samOut.WriteRecord(samHeader, samRecord))
      {
         // Failed to write a record.
         fprintf(stderr, "%s\n", samOut.GetStatusMessage());
         writeStatus = samOut.GetStatus();
      }
   }

   if(samIn.GetStatus() != SamStatus::NO_MORE_RECS)
   {
      // Failed to read a record.
      fprintf(stderr, "%s\n", samIn.GetStatusMessage());
   }

   std::cout << std::endl << "Number of records read = " << 
      samIn.GetCurrentRecordCount() << std::endl;
   std::cout << "Number of records written = " << 
      samOut.GetCurrentRecordCount() << std::endl;

   if(samIn.GetStatus() != SamStatus::NO_MORE_RECS)
   {
      // Failed reading a record.
      return(samIn.GetStatus());
   }

   // Since the reads were successful, return the status based
   // on the status of the writes.  If any failed, return
   // their failure status.
   return(writeStatus);
}

Simple Read/Write Indexed BAM Files

This example reads in the inputFilename bam file and writes it back out section by section to the specified outputFilename, starting with section -1. It also prints a count of the number of records in each section.

#include "SamFile.h"
#include "SamValidation.h"

int ReadIndexedBam(const char* inputFilename,
                   const char* outputFilename,
                   const char* indexFilename)
{
   // Open the input file for reading.
   SamFile samIn;
   if(!samIn.OpenForRead(inputFilename))
   {
      fprintf(stderr, "%s\n", samIn.GetStatusMessage());
      return(samIn.GetStatus());
   }

   // Open the bam index file for reading.
   if(!samIn.ReadBamIndex(indexFilename))
   {
      fprintf(stderr, "%s\n", samIn.GetStatusMessage());
      return(samIn.GetStatus());
   }

   // Open the output file for writing.
   SamFile samOut;
   if(!samOut.OpenForWrite(outputFilename))
   {
      fprintf(stderr, "%s\n", samOut.GetStatusMessage());
      return(samOut.GetStatus());
   }

   // Read the sam header.
   SamFileHeader samHeader;
   if(!samIn.ReadHeader(samHeader))
   {
      fprintf(stderr, "%s\n", samIn.GetStatusMessage());
      return(samIn.GetStatus());
   }

   // Write the sam header.
   if(!samOut.WriteHeader(samHeader))
   {
      fprintf(stderr, "%s\n", samOut.GetStatusMessage());
      return(samOut.GetStatus());     
   }

   SamRecord samRecord;
   SamValidationErrors samInvalidErrors;

   // Loop through each Reference.
   for(int i = -1; i < 23; i++)
   {
      int numSectionRecords = 0;
      samIn.SetReadSection(i);
      // Keep reading records until they aren't more.
      while(samIn.ReadRecord(samHeader, samRecord))
      {
         numSectionRecords++;
         // Successfully read a record from the file, so write it.
         if(!samOut.WriteRecord(samHeader, samRecord))
         {
            // Failed to write a record.
            fprintf(stderr, "%s\n", samOut.GetStatusMessage());
         }
      }

      if(samIn.GetStatus() != SamStatus::NO_MORE_RECS)
      {
         // Failed to read a record.
         fprintf(stderr, "%s\n", samIn.GetStatusMessage());
      }

      std::cout << "Reference ID " << i << " has " << numSectionRecords 
                << " records" << std::endl;
   }
   
   std::cout << "Number of records = " << samIn.GetCurrentRecordCount() << std::endl;
   
   return(0);
}

Validate SAM/BAM Files Continuing Processing After Allowed Errors

This example demonstrates calling validation on a SAM/BAM file. It also demonstrates a more complicated read loop that continues until a configurable number of errors are hit.

#include "SamFile.h"
#include "SamValidation.h"

int Validate(String inFile, int quitAfterErrorNum, int maxReportedErrors, SamFile::SortedType sortType)
{
   // Check to see if the in file was specified.  If it wasn't this is an error.
   if(inFile == "")
   {
      // In file was not specified but it is mandatory.
      std::cerr << "--in is a mandatory argument for validate, "
                << "but was not specified" << std::endl;
      return(-1);
   }

   SamFile samIn;
   // Open the file for reading.   
   if(!samIn.OpenForRead(inFile))
   {
      fprintf(stderr, "%s\n", samIn.GetStatusMessage());
      return(samIn.GetStatus());
   }

   // Set the sorting validation type.
   samIn.setSortedValidation(sortType);

   // Read the sam header.
   SamFileHeader samHeader;
   if(!samIn.ReadHeader(samHeader))
   {
      fprintf(stderr, "%s\n", samIn.GetStatusMessage());
      return(samIn.GetStatus());
   }

   // Read the sam records.
   SamRecord samRecord;

   // Track the status.
   SamStatus::Status status = SamStatus::SUCCESS;

   // Keep reading records until the end of the file is reached.
   int numValidRecords = 0;
   int numInvalidRecords = 0;
   int numErrorRecords = 0;
   int numRecords = 0;
   int numReportedErrors = 0;
   int totalErrorRecords = 0;

   SamValidationErrors invalidSamErrors;

   // Keep reading records from the file until SamFile::ReadRecord
   // indicates to stop (returns false).
   while( (samIn.ReadRecord(samHeader, samRecord))||
          (SamStatus::isContinuableStatus(samIn.GetStatus()) &&
           ( (quitAfterErrorNum < 0) ||
             (numErrorRecords < quitAfterErrorNum) )) )
   {
      ++numRecords;
      if(samIn.GetStatus() == SamStatus::SUCCESS)
      {
         // Successfully set the record, so check to see if it is valid.
         // Clear any errors in the list.
         invalidSamErrors.clear();
         if(!SamValidator::isValid(samHeader, samRecord, invalidSamErrors))
         {
            // The record is not valid.
            ++numInvalidRecords;
            ++totalErrorRecords;
            if(numReportedErrors < maxReportedErrors)
            {
               std::cout << "Record " << numRecords << std::endl
                         << invalidSamErrors << std::endl;
               numReportedErrors += invalidSamErrors.numErrors();
               status = SamStatus::INVALID;
            }
         }
         else
         {
            // Valid record, so increment the counter.
            ++numValidRecords;
         }
      }
      else
      {
         // Error reading the record.
         ++numErrorRecords;
         if(numReportedErrors < maxReportedErrors)
         {
            // report error.
            std::cout << "Record " << numRecords << std::endl
                      << samIn.GetStatusMessage() << std::endl
                      << std::endl;
            ++numReportedErrors;
            ++totalErrorRecords;
            status = samIn.GetStatus();
         }
      }
   }

   if(samIn.GetStatus() != SamStatus::NO_MORE_RECS)
   {
      ++numErrorRecords;
      // The last read call had a failure, so report it.
      if(numReportedErrors < maxReportedErrors)
      {
         std::cout << "Record " << numRecords << ": ";
         std::cout << std::endl << samIn.GetStatusMessage() << std::endl;
      }
      status = samIn.GetStatus();
   }

   printf("\nNumber of records read = %d\n", numRecords);
   printf("Number of valid records = %d\n", numValidRecords);
   printf("Returning: %d (%s)\n", status, SamStatus::getStatusString(status));
   return(status);
}