BamTrimMain.cpp

   1 // ***************************************************************************\r
   2 // BamTrimMain.cpp (c) 2009 Derek Barnett\r
   3 // Marth Lab, Department of Biology, Boston College\r
   4 // All rights reserved.\r
   5 // ---------------------------------------------------------------------------\r
   6 // Last modified: 15 July 2009 (DB)\r
   7 // ---------------------------------------------------------------------------\r
   8 // Basic example of reading/writing BAM files. Pulls alignments overlapping \r
   9 // the range specified by user from one BAM file and writes to a new BAM file.\r
  10 // ***************************************************************************\r
  11 \r
  12 // Std C/C++ includes\r
  13 #include <cstdlib>\r
  14 #include <iostream>\r
  15 #include <string>\r
  16 using namespace std;\r
  17 \r
  18 // BamTools includes\r
  19 #include "BamReader.h"\r
  20 #include "BamWriter.h"\r
  21 using namespace BamTools;\r
  22 \r
  23 int main(int argc, char* argv[]) {\r
  24 \r
  25         // validate argument count\r
  26         if( argc != 7 ) {\r
  27                 cerr << "USAGE: " << argv[0] << " <input BAM file> <input BAM index file> <output BAM file> <reference name> <leftBound> <rightBound> " << endl;\r
  28                 exit(1);\r
  29         }\r
  30 \r
  31         // store arguments\r
  32         string inBamFilename  = argv[1];\r
  33         string indexFilename  = argv[2];\r
  34         string outBamFilename = argv[3];\r
  35         string referenceName  = argv[4];\r
  36         string leftBound_str  = argv[5];\r
  37         string rightBound_str = argv[6];\r
  38 \r
  39         // open our BAM reader\r
  40         BamReader reader;\r
  41         reader.Open(inBamFilename, indexFilename);\r
  42 \r
  43         // get header & reference information\r
  44         string header = reader.GetHeaderText();\r
  45         RefVector references = reader.GetReferenceData();\r
  46         \r
  47         // open our BAM writer\r
  48         BamWriter writer;\r
  49         writer.Open(outBamFilename, header, references);\r
  50 \r
  51         // get reference ID from name\r
  52         int refID = 0;\r
  53         RefVector::const_iterator refIter = references.begin();\r
  54         RefVector::const_iterator refEnd  = references.end();\r
  55         for ( ; refIter != refEnd; ++refIter ) {\r
  56                 if ( (*refIter).RefName == referenceName ) { break; }\r
  57                 ++refID;\r
  58         }\r
  59         \r
  60         // validate ID\r
  61         if ( refIter == refEnd ) {\r
  62                 cerr << "Reference: " << referenceName << " not found." << endl;\r
  63                 exit(1);\r
  64         }\r
  65 \r
  66         // convert boundary arguments to numeric values\r
  67         unsigned int leftBound  = (unsigned int) atoi( leftBound_str.c_str()  );\r
  68         unsigned int rightBound = (unsigned int) atoi( rightBound_str.c_str() );\r
  69         \r
  70         // attempt jump to range of interest\r
  71         if ( reader.Jump(refID, leftBound) ) {\r
  72                 \r
  73                 // while data exists and alignment begin before right bound\r
  74                 BamAlignment bAlignment;\r
  75                 while ( reader.GetNextAlignment(bAlignment) && (bAlignment.Position <= rightBound) ) {\r
  76                         // save alignment to archive\r
  77                         writer.SaveAlignment(bAlignment);\r
  78                 }\r
  79         } \r
  80         \r
  81         // if jump failed\r
  82         else {\r
  83                 cerr << "Could not jump to ref:pos " << referenceName << ":" << leftBound << endl;\r
  84                 exit(1);\r
  85         }\r
  86 \r
  87         // clean up and exit\r
  88         reader.Close();\r
  89         writer.Close(); \r
  90         return 0;\r
  91 }