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
12 // Std C/C++ includes
\r
16 using namespace std;
\r
18 // BamTools includes
\r
19 #include "BamReader.h"
\r
20 #include "BamWriter.h"
\r
21 using namespace BamTools;
\r
23 int main(int argc, char* argv[]) {
\r
25 // validate argument count
\r
27 cerr << "USAGE: " << argv[0] << " <input BAM file> <input BAM index file> <output BAM file> <reference name> <leftBound> <rightBound> " << endl;
\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
39 // open our BAM reader
\r
41 reader.Open(inBamFilename, indexFilename);
\r
43 // get header & reference information
\r
44 string header = reader.GetHeaderText();
\r
45 RefVector references = reader.GetReferenceData();
\r
47 // open our BAM writer
\r
49 writer.Open(outBamFilename, header, references);
\r
51 // get reference ID from name
\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
61 if ( refIter == refEnd ) {
\r
62 cerr << "Reference: " << referenceName << " not found." << endl;
\r
66 // convert boundary arguments to numeric values
\r
67 int leftBound = (int) atoi( leftBound_str.c_str() );
\r
68 int rightBound = (int) atoi( rightBound_str.c_str() );
\r
70 // attempt jump to range of interest
\r
71 if ( reader.Jump(refID, leftBound) ) {
\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
83 cerr << "Could not jump to ref:pos " << referenceName << ":" << leftBound << endl;
\r
87 // clean up and exit
\r