int64_t AlignmentsBeginOffset;\r
string Filename;\r
string IndexFilename;\r
-
+ \r
// system data\r
bool IsBigEndian;\r
\r
int BamReader::GetReferenceCount(void) const { return d->References.size(); }\r
const RefVector BamReader::GetReferenceData(void) const { return d->References; }\r
int BamReader::GetReferenceID(const string& refName) const { return d->GetReferenceID(refName); }\r
-const std::string BamReader::GetFilename(void) const { return d->Filename; }
+const std::string BamReader::GetFilename(void) const { return d->Filename; }\r
\r
// index operations\r
bool BamReader::CreateIndex(void) { return d->CreateIndex(); }\r
const uint64_t& lastOffset)\r
{\r
// get converted offsets\r
- int beginOffset = bAlignment.Position >> BAM_LIDX_SHIFT;
+ int beginOffset = bAlignment.Position >> BAM_LIDX_SHIFT;\r
int endOffset = (bAlignment.GetEndPosition() - 1) >> BAM_LIDX_SHIFT;\r
\r
// resize vector if necessary\r
// read starts after left boundary\r
if ( bAlignment.Position >= CurrentLeft) { return true; }\r
\r
- // return whether alignment end overlaps left boundary
+ // return whether alignment end overlaps left boundary\r
return ( bAlignment.GetEndPosition() >= CurrentLeft );\r
}\r
\r
CXX= g++\r
CXXFLAGS= -Wall -O3\r
-PROG= BamConversion BamDump BamTrim bamtools\r
+PROG= bamtools bamtools_test\r
LIBS= -lz\r
+OBJS= BGZF.o BamReader.o BamWriter.o BamMultiReader.o bamtools.o
\r
all: $(PROG)\r
\r
-bamtools: BGZF.o BamReader.o BamWriter.o BamMultiReader.o bamtools.o\r
- $(CXX) $(CXXFLAGS) -o $@ BGZF.o BamReader.o BamWriter.o BamMultiReader.o bamtools.o $(LIBS)\r
-\r
-BamConversion: BGZF.o BamReader.o BamWriter.o BamConversionMain.o\r
- $(CXX) $(CXXFLAGS) -o $@ BGZF.o BamReader.o BamWriter.o BamConversionMain.o $(LIBS)\r
-\r
-BamDump: BGZF.o BamReader.o BamDumpMain.o\r
- $(CXX) $(CXXFLAGS) -o $@ BGZF.o BamReader.o BamDumpMain.o $(LIBS)\r
-\r
-BamTrim: BGZF.o BamReader.o BamWriter.o BamTrimMain.o\r
- $(CXX) $(CXXFLAGS) -o $@ BGZF.o BamReader.o BamWriter.o BamTrimMain.o $(LIBS)\r
+bamtools: $(OBJS)\r
+ $(CXX) $(CXXFLAGS) -o $@ $(OBJS) $(LIBS)\r
\r
clean:\r
- rm -fr gmon.out *.o *.a a.out *~ $(PROG)\r
+ rm -fr gmon.out *.o *.a a.out *~\r
README : BAMTOOLS
------------------------------------------------------------
-BamTools: a C++ API for reading/writing BAM files.
+BamTools: a C++ API & toolkit for reading/writing/manipulating BAM files.
+
+I. Introduction
+ a. The API
+ b. The Toolkit
+
+II. Usage
+ a. The API
+ b. The Toolkit
-I. Introduction
-II. Usage
III. Contact
------------------------------------------------------------
I. Introduction:
+BamTools provides both a programmer's API and an end-user's toolkit for handling
+BAM files.
+Ia. The API
The API consists of 2 main modules - BamReader and BamWriter. As you would expect,
BamReader provides read-access to BAM files, while BamWriter does the writing of BAM
files. BamReader provides an interface for random-access (jumping) in a BAM file,
as well as generating BAM index files.
+BamMultiReader is an extra module that allows you to manage multiple open BAM file
+for reading. It provides some validation & bookkeeping under the hood to keep all
+files sync'ed for
+
An additional file, BamAux.h, is included as well.
This file contains the common data structures and typedefs used throught the API.
BGZF.h & BGZF.cpp contain our implementation of the Broad Institute's
BGZF compression format.
-BamConversion, BamDump, and BamTrim are 3 'toy' examples on how the API can be used.
+
+Ib. The Toolkit
+If you've been using BamTools since the early days, you'll notice that our 'toy' API
+examples (BamConversion, BamDump, and BamTrim) are now gone. In their place is a set
+of features we hope you find useful.
+
+** More explanation here **
+
+usage: bamtools [--help] COMMAND [ARGS]
+
+Available bamtools commands:
+ coverage Prints coverage statistics from the input BAM file
+ dump Dump BAM file contents to text output
+ header Prints BAM header information
+ index Generates index for BAM file
+ merge Merge multiple BAM files into single file
+ sort Sorts the BAM file according to some criteria
+ stats Prints some basic statistics from the input BAM file
+
+See 'bamtools help COMMAND' for more information on a specific command.
+
+** Follow-up explanation here **
------------------------------------------------------------
II. Usage :
+** General usage information - perhaps explain common terms, point to SAM/BAM spec, etc **
+
+
+IIa. The API
+
To use this API, you simply need to do 3 things:
1 - Drop the BamTools files somewhere the compiler can find them.
See any included programs and Makefile for more specific compiling/usage examples.
See comments in the header files for more detailed API documentation.
+
+IIb. The Toolkit
+
+** More indepth overview for the toolkit commands **
+
------------------------------------------------------------
III. Contact :
Feel free to contact me with any questions, comments, suggestions, bug reports, etc.
- Derek Barnett
+Marth Lab
+Biology Dept., Boston College
+
Email: barnetde@bc.edu
-Project Website: http://sourceforge.net/projects/bamtools
+Project Websites: http://github.com/pezmaster31/bamtools (ACTIVE SUPPORT)
+ http://sourceforge.net/projects/bamtools (major updates only)
// ***************************************************************************
-// bamtools.cpp (c) 2010 Erik Garrison
+// bamtools.cpp (c) 2010 Derek Barnett, Erik Garrison
// Marth Lab, Department of Biology, Boston College
// All rights reserved.
// ---------------------------------------------------------------------------
+// Last modified: 26 May 2010
+// ---------------------------------------------------------------------------
// Integrates a number of BamTools functionalities into a single executable.
// ***************************************************************************
// Std C/C++ includes
-#include <cstdlib>
#include <iostream>
-#include <fstream>
-#include <string>
-using namespace std;
// BamTools includes
-#include "BamReader.h"
-#include "BamWriter.h"
-#include "BamMultiReader.h"
-using namespace BamTools;
-
-void usageSummary(string cmdname) {
- cerr << "usage: " << cmdname << " <command> [options]" << endl
- << "actions:" << endl
- << " index <BAM file> # generates BAM index <BAM file>.bai" << endl
- << " merge <merged BAM file> [<BAM file> <BAM file> ...] # merges BAM files into a single file" << endl
- << " dump [<BAM file> <BAM file> ...] # dumps alignment summaries to stdout" << endl
- << " header [<BAM file> <BAM file> ...] # prints header, or unified header for BAM file or files" << endl;
- //<< " trim <input BAM file> <input BAM index file> <output BAM file> <reference name> <leftBound> <rightBound>" << endl;
-}
-
+#include "bamtools_coverage.h"
+#include "bamtools_dump.h"
+#include "bamtools_header.h"
+#include "bamtools_index.h"
+#include "bamtools_merge.h"
+#include "bamtools_sam.h"
+#include "bamtools_sort.h"
+#include "bamtools_stats.h"
-void BamMerge(string outputFilename, vector<string> filenames) {
-
- BamMultiReader reader;
-
- reader.Open(filenames, false); // opens the files without checking for indexes
-
- string mergedHeader = reader.GetHeaderText();
-
- RefVector references = reader.GetReferenceData();
-
- // open BamWriter
- BamWriter writer;
- writer.Open( outputFilename.c_str(), mergedHeader, references);
+using namespace std;
+using namespace BamTools;
- BamAlignment bAlignment;
- while (reader.GetNextAlignment(bAlignment)) {
- // write to output file
- writer.SaveAlignment(bAlignment);
+// ------------------------------------------
+// bamtools subtool names
+static const string COVERAGE = "coverage";
+static const string DUMP = "dump"; // <-- do we even want to keep this? I think 'bamtools sam' will be more useful anyway
+ // nobody's going to want what was essentially an early, bloated, debugging output
+static const string HEADER = "header";
+static const string INDEX = "index";
+static const string MERGE = "merge";
+static const string SAM = "sam";
+static const string SORT = "sort";
+static const string STATS = "stats";
+
+// ------------------------------------------
+// bamtools help/version names
+static const string HELP = "help";
+static const string LONG_HELP = "--help";
+static const string SHORT_HELP = "-h";
+
+static const string VERSION = "version";
+static const string LONG_VERSION = "--version";
+static const string SHORT_VERSION = "-v";
+
+int Help(int argc, char* argv[]) {
+
+ // 'bamtools help COMMAND'
+ if (argc > 2) {
+ if ( argv[2] == COVERAGE) return BamCoverageHelp();
+ if ( argv[2] == DUMP ) return BamDumpHelp(); // keep?
+ if ( argv[2] == HEADER ) return BamHeaderHelp();
+ if ( argv[2] == INDEX ) return BamIndexHelp();
+ if ( argv[2] == MERGE ) return BamMergeHelp();
+ if ( argv[2] == SAM ) return BamSamHelp();
+ if ( argv[2] == SORT ) return BamSortHelp();
+ if ( argv[2] == STATS ) return BamStatsHelp();
}
-
- // close output file
- writer.Close();
- // close input files
- reader.Close();
-
-}
-
-void BamCreateIndex(const char* inputFilename) {
-
- // open our BAM reader
- BamReader reader;
- reader.Open(inputFilename);
-
- // create index file
- reader.CreateIndex();
-
- // close our file
- reader.Close();
-
-}
-
-// Spit out basic BamAlignment data
-void PrintAlignment(const BamAlignment& alignment) {
- cout << "---------------------------------" << endl;
- cout << "Name: " << alignment.Name << endl;
- cout << "Aligned to: " << alignment.RefID;
- cout << ":" << alignment.Position << endl;
- cout << endl;
-}
-
-void BamDump(vector<string> files) {
-
- BamMultiReader reader;
- reader.Open(files);
-
- BamAlignment bAlignment;
- while (reader.GetNextAlignment(bAlignment)) {
- PrintAlignment(bAlignment);
- }
-
- reader.Close();
-
+
+ // either 'bamtools help' or unrecognized argument after 'help'
+ cerr << endl;
+ cerr << "usage: bamtools [--help] COMMAND [ARGS]" << endl;
+ cerr << endl;
+ cerr << "Available bamtools commands:" << endl;
+ cerr << "\tcoverage\tPrints coverage statistics from the input BAM file" << endl;
+ cerr << "\tdump\t\tDump BAM file contents to text output" << endl; // keep?
+ cerr << "\theader\t\tPrints BAM header information" << endl;
+ cerr << "\tindex\t\tGenerates index for BAM file" << endl;
+ cerr << "\tmerge\t\tMerge multiple BAM files into single file" << endl;
+ cerr << "\tsam\t\tPrints the BAM file in SAM (text) format" << endl;
+ cerr << "\tsort\t\tSorts the BAM file according to some criteria" << endl;
+ cerr << "\tstats\t\tPrints some basic statistics from the input BAM file" << endl;
+ cerr << endl;
+ cerr << "See 'bamtools help COMMAND' for more information on a specific command." << endl;
+ cerr << endl;
+
+ return 0;
}
-void BamDumpHeader(vector<string> files) {
-
- BamMultiReader reader;
- reader.Open(files, false);
-
- cout << reader.GetHeaderText() << endl;
-
- reader.Close();
-
+int Version(void) {
+ cout << endl;
+ cout << "bamtools v0.x.xx" << endl;
+ cout << "Part of BamTools API and toolkit" << endl;
+ cout << "Primary authors: Derek Barnett, Erik Garrison, Michael Stromberg" << endl;
+ cout << "(c) 2009-2010 Marth Lab, Biology Dept., Boston College" << endl;
+ cout << endl;
+ return 0;
}
-
int main(int argc, char* argv[]) {
- // validate argument count
- if( argc < 3 ) {
- usageSummary(argv[0]);
- exit(1);
- }
-
- string command = argv[1];
+ // just 'bamtools'
+ if ( (argc == 1) ) return Help(argc, argv);
- if (command == "index") {
-
- BamCreateIndex(argv[2]);
-
- } else if (command == "merge") {
-
- vector<string> files;
- string outputFile = argv[2];
-
- // check if our output exists, and exit if so
- ifstream output(outputFile.c_str());
- if (output.good()) {
- cerr << "ERROR: output file " << outputFile << " exists, exiting." << endl;
- exit(1);
- } else {
- output.close();
- }
-
- for (int i = 3; i<argc; ++i) {
- files.push_back(argv[i]);
- }
- BamMerge(outputFile, files);
-
- } else if (command == "dump") {
+ // 'bamtools help', 'bamtools --help', or 'bamtools -h'
+ if ( (argv[1] == HELP) || (argv[1] == LONG_HELP) || (argv[1] == SHORT_HELP) ) return Help(argc, argv);
+
+ // 'bamtools version', 'bamtools --version', or 'bamtools -v'
+ if ( (argv[1] == VERSION) || (argv[1] == LONG_VERSION) || (argv[1] == SHORT_VERSION) ) return Version();
- vector<string> files;
- for (int i = 2; i<argc; ++i) {
- files.push_back(argv[i]);
- }
- BamDump(files);
-
- } else if (command == "header") {
-
- vector<string> files;
- for (int i = 2; i<argc; ++i) {
- files.push_back(argv[i]);
- }
- BamDumpHeader(files);
-
- }
-
-
- return 0;
+ // run desired sub-tool
+ if ( argv[1] == COVERAGE ) return RunBamCoverage(argc, argv);
+ if ( argv[1] == DUMP ) return RunBamDump(argc, argv); // keep?
+ if ( argv[1] == HEADER ) return RunBamHeader(argc, argv);
+ if ( argv[1] == INDEX ) return RunBamIndex(argc, argv);
+ if ( argv[1] == MERGE ) return RunBamMerge(argc, argv);
+ if ( argv[1] == SAM ) return RunBamSam(argc, argv);
+ if ( argv[1] == SORT ) return RunBamSort(argc, argv);
+ if ( argv[1] == STATS ) return RunBamStats(argc, argv);
+
+ // unrecognized 2nd argument, print help
+ return Help(argc, argv);
}
--- /dev/null
+// ***************************************************************************
+// bamtools_coverage.h (c) 2010 Derek Barnett, Erik Garrison
+// Marth Lab, Department of Biology, Boston College
+// All rights reserved.
+// ---------------------------------------------------------------------------
+// Last modified: 26 May 2010
+// ---------------------------------------------------------------------------
+// Prints coverage statistics for a single BAM file
+//
+// ** Expand to multiple??
+//
+// ***************************************************************************
+
+#ifndef BAMTOOLS_COVERAGE_H
+#define BAMTOOLS_COVERAGE_H
+
+#include <iostream>
+#include <string>
+
+#include "BamReader.h"
+#include "bamtools_getopt.h"
+
+namespace BamTools {
+
+int BamCoverageHelp(void) {
+ std::cerr << std::endl;
+ std::cerr << "usage:\tbamtools coverage [--in BAM file]" << std::endl;
+ std::cerr << "\t-i, --in\tInput BAM file to generate coverage stats\t[default=stdin]" << std::endl;
+ std::cerr << std::endl;
+ return 0;
+}
+
+int RunBamCoverage(int argc, char* argv[]) {
+
+ // else parse command line for args
+ GetOpt options(argc, argv, 1);
+
+ std::string inputFilename;
+ options.addOption('i', "in", &inputFilename);
+
+ if ( !options.parse() ) return BamCoverageHelp();
+ if ( inputFilename.empty() ) { inputFilename = "stdin"; }
+
+// // open our BAM reader
+// BamReader reader;
+// reader.Open(inputFilename);
+
+ // generate coverage stats
+ std::cerr << "Generating coverage stats for " << inputFilename << std::endl;
+ std::cerr << "FEATURE NOT YET IMPLEMENTED!" << std::endl;
+
+ // clean & exit
+// reader.Close();
+ return 0;
+}
+
+} // namespace BamTools
+
+#endif // BAMTOOLS_COVERAGE_H
--- /dev/null
+// ***************************************************************************
+// bamtools_dump.h (c) 2010 Derek Barnett, Erik Garrison
+// Marth Lab, Department of Biology, Boston College
+// All rights reserved.
+// ---------------------------------------------------------------------------
+// Last modified: 26 May 2010
+// ---------------------------------------------------------------------------
+// Dumps alignment summaries out to stdout.
+//
+// ** This should probably go the way of the dodo soon? bamtools sam makes this
+// obsolete and probably worthless.
+//
+// ***************************************************************************
+
+#ifndef BAMTOOLS_DUMP_H
+#define BAMTOOLS_DUMP_H
+
+#include <iostream>
+#include <string>
+#include <vector>
+
+#include "BamMultiReader.h"
+// #include "GetOpt.h"
+#include "bamtools_getopt.h"
+
+namespace BamTools {
+
+int BamDumpHelp(void) {
+ std::cerr << std::endl;
+ std::cerr << "usage:\tbamtools dump [BAM file1] [BAM file2] [BAM file3]..." << std::endl;
+ std::cerr << "\t[BAM file]\tInput file(s) to dump alignment summaries from [default=stdin]" << std::endl;
+ std::cerr << std::endl;
+ return 0;
+}
+
+// Spit out basic BamAlignment data
+void PrintAlignment(const BamTools::BamAlignment& alignment) {
+ std::cout << "---------------------------------" << std::endl;
+ std::cout << "Name: " << alignment.Name << std::endl;
+ std::cout << "Aligned to: " << alignment.RefID;
+ std::cout << ":" << alignment.Position << std::endl;
+ std::cout << std::endl;
+}
+
+int RunBamDump(int argc, char* argv[]) {
+
+ // else parse command line for args
+ GetOpt options(argc, argv, 1);
+
+ std::vector<std::string> inputFilenames;
+ options.addVariableLengthOption("in", &inputFilenames);
+
+ if ( !options.parse() ) return BamDumpHelp();
+ if ( inputFilenames.empty() ) { inputFilenames.push_back("stdin"); }
+
+ // open files
+ BamMultiReader reader;
+ reader.Open(inputFilenames, false);
+
+ // dump alignment summaries to stdout
+ BamAlignment bAlignment;
+ while (reader.GetNextAlignment(bAlignment)) {
+ PrintAlignment(bAlignment);
+ }
+
+ // clean up & exit
+ reader.Close();
+ return 0;
+}
+
+} // namespace BamTools
+
+#endif // BAMTOOLS_DUMP_H
--- /dev/null
+// ***************************************************************************
+// bamtools_getopt.h (c) 2010 Derek Barnett, Erik Garrison
+// Marth Lab, Department of Biology, Boston College
+// All rights reserved.
+// ---------------------------------------------------------------------------
+// Last modified: 26 May 2010
+// ---------------------------------------------------------------------------
+// Provides a configurable commandline parser used by the BamTools subtools
+// ***************************************************************************
+
+#ifndef BAMTOOLS_GETOPT_H
+#define BAMTOOLS_GETOPT_H
+
+// C includes
+#include <cassert>
+#include <cstdlib>
+
+// C++ includes
+#include <iostream>
+#include <map>
+#include <string>
+#include <vector>
+
+namespace BamTools {
+
+class GetOpt {
+
+ // ctors & dtor
+ public:
+
+ // ctor: takes the 'standard' command line args (optional offset)
+ GetOpt(int argc, char* argv[], int offset = 0);
+
+ // d-tor
+ ~GetOpt(void);
+
+ // set rules for bare word arguments
+ public:
+ // add an optional 'bare word' argument (eg 'help')
+ // 'name' is not used on the command line, but for reporting
+ void addOptionalArgument(const std::string& name, std::string* value);
+
+ // add a required 'bare word' argument (eg input data file)
+ // 'name' is not used on the command line, but for reporting
+ void addRequiredArgument(const std::string& name, std::string* value);
+
+ // set rules for key=>value options
+ public:
+ // add standard option with arguments ( -Wall, -O2, --type=foo )
+ void addOption(const char shortName, const std::string& longName, std::string* value);
+
+ // add an option whose argument is optional (eg --log may default to dumping to stderr, unless a file is specified )
+ // must provide a default string
+ void addOptionalOption(const char shortName, const std::string& longName, std::string* value, const std::string& defaultValue);
+ void addOptionalOption(const std::string& longName, std::string* value, const std::string& defaultValue);
+
+ // add a repeatable option (like compiler includes -I/path/ -I/path2/ etc)
+ // only supporting one type of name (short/long) for this option for now
+ void addRepeatableOption(const char shortName, std::vector<std::string>* values); // single char version
+ void addRepeatableOption(const std::string& longName, std::vector<std::string>* values); // long name version
+
+ // add an option that takes a variable number of arguments ( --files f1 f2 f3 f4... )
+ void addVariableLengthOption(const std::string& longName, std::vector<std::string>* values);
+
+ // set rules for on/off switch
+ public:
+ // on/off switch ( --verbose --searchOnly ) only long names supported for now
+ void addSwitch(const std::string& longName, bool* ok);
+
+ // parse and query methods
+ public:
+
+ // get application name
+ const std::string& applicationName(void) const;
+
+ // query if particular 'bare-word' argument is set
+ bool isSet(const std::string& name) const;
+
+ // runs parser (does validation and assign values to arguments)
+ // returns success/fail
+ bool parse(void);
+
+ void print(void);
+
+ // define Option-related types & enums
+ private:
+ enum OptionType { OptUnknown = 0
+ , OptEnd
+ , OptSwitch
+ , OptArg1
+ , OptOptional
+ , OptRepeat
+ , OptVariable
+ };
+
+ // define Option
+ struct Option {
+
+ // ctor
+ Option(OptionType t = OptUnknown, const char shortName = 0, const std::string& longName = "")
+ : Type(t)
+ , ShortName(shortName)
+ , LongName(longName)
+ , BoolValue(0)
+ { }
+
+ // data members
+ OptionType Type;
+ char ShortName;
+ std::string LongName;
+ union {
+ bool* BoolValue;
+ std::string* StringValue;
+ std::vector<std::string>* ListValue;
+ };
+ std::string Default;
+ };
+
+ // internal methods
+ private:
+ void init(int argc, char* argv[], int offset);
+ void saveOption(const Option& opt); // const & ?? he doesnt use it - why?
+ void setSwitch(const Option& opt);
+
+ // data members
+ private:
+ std::vector<Option> m_options;
+ std::map<std::string, int> m_setOptions;
+ std::vector<std::string> m_args;
+ std::string m_appname;
+
+ int m_numberRequiredArguments;
+ int m_numberOptionalArguments;
+ Option m_requiredArgument;
+ Option m_optionalArgument;
+
+ int m_currentArgument;
+};
+
+inline
+GetOpt::GetOpt(int argc, char* argv[], int offset)
+{
+ init(argc, argv, offset);
+}
+
+inline
+GetOpt::~GetOpt(void) { }
+
+// add an optional 'bare word' argument (eg 'help')
+// 'name' is not used on the command line, but for reporting
+inline
+void GetOpt::addOptionalArgument(const std::string& name, std::string* value) {
+
+ Option opt( OptUnknown, 0, name );
+ opt.StringValue = value;
+ m_optionalArgument = opt;
+ ++m_numberOptionalArguments;
+ *value = std::string();
+}
+
+// add a required 'bare word' argument (eg input data file)
+// 'name' is not used on the command line, but for reporting
+inline
+void GetOpt::addRequiredArgument(const std::string& name, std::string* value) {
+
+ Option opt( OptUnknown, 0, name );
+ opt.StringValue = value;
+ m_requiredArgument = opt;
+ ++m_numberRequiredArguments;
+ *value = std::string();
+}
+
+// add standard option with arguments ( -Wall, -O2, --type=foo )
+inline
+void GetOpt::addOption(const char shortName, const std::string& longName, std::string* value) {
+
+ Option opt( OptArg1, shortName, longName );
+ opt.StringValue = value;
+ saveOption(opt);
+ *value = std::string();
+}
+
+// add an option whose argument is optional (eg --log may default to dumping to stderr, unless a file is specified )
+// must provide a default string
+// short & long name version
+inline
+void GetOpt::addOptionalOption(const char shortName, const std::string& longName, std::string* value, const std::string& defaultValue) {
+
+ Option opt( OptOptional, shortName, longName );
+ opt.StringValue = value;
+ opt.Default = defaultValue;
+ saveOption(opt);
+ *value = std::string();
+}
+
+// long name only version
+inline
+void GetOpt::addOptionalOption(const std::string& longName, std::string* value, const std::string& defaultValue) {
+ addOptionalOption(0, longName, value, defaultValue);
+}
+
+// add a repeatable option (like compiler includes -I/path/ -I/path2/ etc)
+// only supporting one type of name (short/long) for this option for now
+// short name only version
+inline
+void GetOpt::addRepeatableOption(const char shortName, std::vector<std::string>* values) {
+
+ Option opt( OptRepeat, shortName, std::string() );
+ opt.ListValue = values;
+ saveOption(opt);
+ *values = std::vector<std::string>();
+}
+
+// long name only version
+inline
+void GetOpt::addRepeatableOption(const std::string& longName, std::vector<std::string>* values) {
+
+ Option opt( OptRepeat, 0, longName );
+ opt.ListValue = values;
+ saveOption(opt);
+ *values = std::vector<std::string>();
+}
+
+// add an option that takes a variable number of arguments ( --files f1 f2 f3 f4... )
+inline
+void GetOpt::addVariableLengthOption(const std::string& longName, std::vector<std::string>* values) {
+
+ Option opt( OptVariable, 0, longName );
+ opt.ListValue = values;
+ saveOption(opt);
+ *values = std::vector<std::string>();
+}
+
+// on/off switch ( --verbose --searchOnly ) only long names supported for now
+inline
+void GetOpt::addSwitch(const std::string& longName, bool* ok) {
+
+ Option opt( OptSwitch, 0, longName );
+ opt.BoolValue = ok;
+ saveOption(opt);
+ *ok = false;
+}
+
+inline
+const std::string& GetOpt::applicationName(void) const {
+ return m_appname;
+}
+
+inline
+void GetOpt::init(int argc, char* argv[], int offset) {
+
+ m_numberRequiredArguments = 0;
+ m_numberOptionalArguments = 0;
+ m_currentArgument = 1;
+
+ if ( argc > 0 ) {
+
+ // store app name
+ std::string fullPath = argv[0];
+ size_t lastSlash = fullPath.find_last_of("/\\"); // should work on Unix- and Windows-style paths
+ m_appname = fullPath.substr(lastSlash + 1);
+
+ // store remaining arguments from offset to end
+ for (int i = offset + 1; i < argc; ++i) {
+ m_args.push_back( argv[i] );
+ }
+
+ } else {
+ std::cerr << "GetOpt ERROR: No arguments given." << std::endl;
+ exit(1);
+ }
+}
+
+// query if particular 'bare-word' argument is set
+inline
+bool GetOpt::isSet(const std::string& name) const {
+ return ( m_setOptions.find(name) != m_setOptions.end() );
+}
+
+// runs parser (does validation and assign values to arguments)
+// returns success/fail
+inline
+bool GetOpt::parse(void) {
+
+ // initialize argument stack (reversed input args)
+ std::vector<std::string> argStack( m_args.rbegin(), m_args.rend() );
+
+ // initialize state
+ enum State { StartingState, ExpectingState, OptionalState };
+ State state = StartingState;
+
+ // initialize token types
+ enum TokenType { LongOpt, ShortOpt, Arg, End };
+ TokenType token = End;
+ TokenType currentType = End;
+
+ // store option list bounds
+ std::vector<Option>::const_iterator optBegin = m_options.begin();
+ std::vector<Option>::const_iterator optEnd = m_options.end();
+
+ // declare currentOption
+ Option currentOption;
+
+ // we're going to fake an 'End' argument
+ bool isExtraLoopNeeded = true;
+
+ // iterate through stack contents & do one extra loop for the fake 'End'
+ while ( !argStack.empty() || isExtraLoopNeeded ) {
+
+ std::string arg;
+ std::string originalArg; // store the original arg because we're going to mangle 'arg'
+
+ // if contents on the arg stack
+ if ( !argStack.empty() ) {
+
+ arg = argStack.back();
+ argStack.pop_back();
+ ++m_currentArgument;
+ originalArg = arg;
+
+ // long option version
+ if ( arg.substr(0,2) == "--" ) {
+
+ // set token type
+ token = LongOpt;
+
+ // strip the '--'
+ arg = arg.substr(2);
+
+ // make sure there's still somthing there
+ if ( arg.empty() ) {
+ std::cerr << "'--' feature is not supported, yet." << std::endl;
+ exit(1);
+ }
+
+ // split any key=value style args
+ size_t foundEqual = arg.find('=');
+ if ( foundEqual != std::string::npos ) {
+
+ // push value back onto stack
+ argStack.push_back( arg.substr(foundEqual+1) );
+ --m_currentArgument;
+
+ // save key as current arg
+ arg = arg.substr(0, foundEqual);
+ }
+
+ }
+
+ // short option version
+ else if ( arg.at(0) == '-' ) {
+
+ // set token type
+ token = ShortOpt;
+
+ // if option is directly followed by argument (eg -Wall), push that arg back onto stack
+ if ( arg.length() > 2 ) {
+ argStack.push_back( arg.substr(2) );
+ --m_currentArgument;
+ }
+
+ // strip the '-'
+ arg = arg[1];
+ }
+
+ // bare-word argument
+ else { token = Arg; }
+ }
+
+ // in fake End iteration
+ else { token = End; }
+
+ // look up arg in list of known options, modify token type if necessary
+ Option opt;
+ if ( token != End ) {
+
+ // look up arg in option list
+ std::vector<Option>::const_iterator optIter = optBegin;
+ for ( ; optIter != optEnd; ++optIter ) {
+ const Option& o = (*optIter);
+ if ( (token == LongOpt && arg == o.LongName) ||
+ (token == ShortOpt && arg.at(0) == o.ShortName) ) {
+ opt = o;
+ break;
+ }
+ }
+
+ // modify token type if needed
+ if ( token == LongOpt && opt.Type == OptUnknown ) {
+ if ( currentOption.Type != OptVariable ) {
+ std::cerr << "GetOpt ERROR: Unknown option --" << arg << std::endl;
+ return false;
+ } else {
+ token = Arg;
+ }
+ } else if ( token == ShortOpt && opt.Type == OptUnknown ) {
+ if ( currentOption.Type != OptVariable ) {
+ std::cerr << "GetOpt ERROR: Unknown option -" << arg.at(0) << std::endl;
+ return false;
+ } else {
+ token = Arg;
+ }
+ }
+ } else { opt = Option(OptEnd); }
+
+
+ // interpret result
+ switch ( state ) {
+
+ case ( StartingState ) :
+
+ if ( opt.Type == OptSwitch ) {
+ setSwitch(opt);
+ m_setOptions.insert( std::pair<std::string, int>(opt.LongName, 1) );
+ m_setOptions.insert( std::pair<std::string, int>((const char*)&opt.ShortName, 1) );
+ } else if ( opt.Type == OptArg1 || opt.Type == OptRepeat ) {
+ state = ExpectingState;
+ currentOption = opt;
+ currentType = token;
+ m_setOptions.insert( std::pair<std::string, int>(opt.LongName, 1) );
+ m_setOptions.insert( std::pair<std::string, int>((const char*)&opt.ShortName, 1) );
+ } else if ( opt.Type == OptOptional || opt.Type == OptVariable ) {
+ state = OptionalState;
+ currentOption = opt;
+ currentType = token;
+ m_setOptions.insert( std::pair<std::string, int>(opt.LongName, 1) );
+ m_setOptions.insert( std::pair<std::string, int>((const char*)&opt.ShortName, 1) );
+ } else if ( opt.Type == OptEnd ) {
+ // do nothing (we're almost done here)
+ } else if ( opt.Type == OptUnknown && token == Arg ) {
+ if ( m_numberRequiredArguments > 0 ) {
+ if ( (*m_requiredArgument.StringValue).empty() ) {
+ *m_requiredArgument.StringValue = arg;
+ } else {
+ std::cerr << "Too many bare arguments" << std::endl;
+ return false;
+ }
+ }
+
+ else if ( m_numberOptionalArguments > 0 ) {
+ if ( (*m_optionalArgument.StringValue).empty() ) {
+ *m_optionalArgument.StringValue = arg;
+ } else {
+ std::cerr << "Too many bare arguments" << std::endl;
+ return false;
+ }
+ }
+ } else {
+ std::cerr << "GetOpt ERROR: Unhandled StartingState case: " << opt.Type << std::endl;
+ exit(1);
+ }
+
+ break;
+
+ case ( ExpectingState ) :
+
+ if ( token == Arg ) {
+ if ( currentOption.Type == OptArg1 ) {
+ *currentOption.StringValue = arg;
+ state = StartingState;
+ } else if ( currentOption.Type == OptRepeat ) {
+ currentOption.ListValue->push_back(arg);
+ state = StartingState;
+ } else {
+ std::cerr << "GetOpt ERROR: Unhandled ExpectingState case: " << currentOption.Type << std::endl;
+ exit(1);
+ }
+ } else {
+ std::string name = (currentType == LongOpt) ? currentOption.LongName : (const char*)¤tOption.ShortName;
+ std::cerr << "GetOpt ERROR: Expected an argument after option: " << name << std::endl;
+ exit(1);
+ }
+
+ break;
+
+ case ( OptionalState ) :
+
+ if ( token == Arg ) {
+ if ( currentOption.Type == OptOptional ) {
+ *currentOption.StringValue = arg;
+ state = StartingState;
+ } else if ( currentOption.Type == OptVariable ) {
+ currentOption.ListValue->push_back(originalArg);
+ // stay in this state
+ } else {
+ std::cerr << "GetOpt ERROR: Unhandled OptionalState case: " << currentOption.Type << std::endl;
+ exit(1);
+ }
+ } else {
+
+ // optional argument not specified
+ if ( currentOption.Type == OptOptional ) {
+ *currentOption.StringValue = currentOption.Default;
+ }
+
+ if ( token != End ) {
+ // re-evaluate current argument
+ argStack.push_back( originalArg );
+ --m_currentArgument;
+ }
+
+ state = StartingState;
+ }
+
+ break;
+ }
+
+ if ( token == End ) {
+ isExtraLoopNeeded = false;
+ }
+ }
+
+ // check that required argument has been satisfied
+ if ( m_numberRequiredArguments > 0 && (*m_requiredArgument.StringValue).empty() ) {
+ std::cerr << "Lacking required argument" << std::endl;
+ return false;
+ }
+
+ return true;
+}
+
+inline
+void GetOpt::print(void) {
+
+ std::cout << "---------------------------------" << std::endl;
+ std::cout << "Options for app: " << m_appname << std::endl;
+ std::cout << std::endl;
+ std::cout << "Args: ";
+ std::vector<std::string>::const_iterator argIter = m_args.begin();
+ std::vector<std::string>::const_iterator argEnd = m_args.end();
+ for ( ; argIter != argEnd; ++argIter ) {
+ std::cout << (*argIter) << " ";
+ }
+ std::cout << std::endl;
+}
+
+inline
+void GetOpt::saveOption(const Option& opt) {
+ // check for conflicts (duplicating options) ??
+ m_options.push_back(opt);
+}
+
+inline
+void GetOpt::setSwitch(const Option& opt) {
+ assert( opt.Type == OptSwitch );
+ *opt.BoolValue = true;
+}
+
+} // namespace BamTools
+
+#endif // BAMTOOLS_GETOPT_H
\ No newline at end of file
--- /dev/null
+// ***************************************************************************
+// bamtools_header.h (c) 2010 Derek Barnett, Erik Garrison
+// Marth Lab, Department of Biology, Boston College
+// All rights reserved.
+// ---------------------------------------------------------------------------
+// Last modified: 26 May 2010
+// ---------------------------------------------------------------------------
+// Prints the SAM-style header from a single BAM file (or merged header from
+// multiple BAM files) to stdout.
+// ***************************************************************************
+
+#ifndef BAMTOOLS_HEADER_H
+#define BAMTOOLS_HEADER_H
+
+#include <iostream>
+#include <string>
+#include <vector>
+
+#include "BamReader.h"
+#include "BamMultiReader.h"
+// #include "GetOpt.h"
+#include "bamtools_getopt.h"
+
+namespace BamTools {
+
+int BamHeaderHelp(void) {
+ std::cerr << std::endl;
+ std::cerr << "usage:\tbamtools header [BAM file1] [BAM file2] [BAM file3]..." << std::endl;
+ std::cerr << "\t[BAM file]\tInput file(s) to dump header contents from [default=stdin]" << std::endl;
+ std::cerr << std::endl;
+ return 0;
+}
+
+int RunBamHeader(int argc, char* argv[]) {
+
+ // else parse command line for args
+ GetOpt options(argc, argv, 1);
+
+ std::vector<std::string> inputFilenames;
+ options.addVariableLengthOption("in", &inputFilenames);
+
+ if ( !options.parse() ) return BamHeaderHelp();
+ if ( inputFilenames.empty() ) { inputFilenames.push_back("stdin"); }
+
+ // open files
+ BamMultiReader reader;
+ reader.Open(inputFilenames, false);
+
+ // dump header contents to stdout
+ std::cout << reader.GetHeaderText() << std::endl;
+
+ // clean up & exit
+ reader.Close();
+ return 0;
+}
+
+} // namespace BamTools
+
+#endif // BAMTOOLS_HEADER_H
--- /dev/null
+// ***************************************************************************
+// bamtools_index.h (c) 2010 Derek Barnett, Erik Garrison
+// Marth Lab, Department of Biology, Boston College
+// All rights reserved.
+// ---------------------------------------------------------------------------
+// Last modified: 26 May 2010
+// ---------------------------------------------------------------------------
+// Creates a BAM index (".bai") file for the provided BAM file.
+// ***************************************************************************
+
+#ifndef BAMTOOLS_INDEX_H
+#define BAMTOOLS_INDEX_H
+
+#include <iostream>
+#include <string>
+
+#include "BamReader.h"
+// #include "GetOpt.h"
+#include "bamtools_getopt.h"
+
+namespace BamTools {
+
+int BamIndexHelp(void) {
+ std::cerr << std::endl;
+ std::cerr << "usage:\tbamtools index [--nclist] <BAM file>" << std::endl;
+ std::cerr << "\t--nclist\tUse NCList indexing scheme (faster?)\t[default=off] ** JUST HERE AS POSSIBLE SWITCH EXAMPLE FOR NOW **" << std::endl;
+ std::cerr << "\t<BAM file>\tInput BAM file to generate index from\t[req'd]" << std::endl;
+ std::cerr << std::endl;
+ return 0;
+}
+
+int RunBamIndex(int argc, char* argv[]) {
+
+ // else parse command line for args
+ GetOpt options(argc, argv, 1);
+
+ std::string inputFilename;
+ options.addRequiredArgument("input", &inputFilename);
+
+ bool useNCList;
+ options.addSwitch("nclist", &useNCList);
+
+ if ( !options.parse() ) return BamIndexHelp();
+
+ // open our BAM reader
+ BamReader reader;
+ reader.Open(inputFilename);
+
+ // create index for BAM file
+ reader.CreateIndex();
+
+ // clean & exit
+ reader.Close();
+ return 0;
+}
+
+} // namespace BamTools
+
+#endif // BAMTOOLS_INDEX_H
--- /dev/null
+// ***************************************************************************
+// bamtools_merge.h (c) 2010 Derek Barnett, Erik Garrison
+// Marth Lab, Department of Biology, Boston College
+// All rights reserved.
+// ---------------------------------------------------------------------------
+// Last modified: 26 May 2010
+// ---------------------------------------------------------------------------
+// Merges multiple BAM files into one.
+//
+// ** Provide selectable region? eg chr2:10000..20000
+//
+// ***************************************************************************
+
+#ifndef BAMTOOLS_MERGE_H
+#define BAMTOOLS_MERGE_H
+
+#include <iostream>
+#include <string>
+#include <vector>
+
+#include "BamMultiReader.h"
+#include "BamWriter.h"
+// #include "GetOpt.h"
+#include "bamtools_getopt.h"
+
+namespace BamTools {
+
+int BamMergeHelp(void) {
+ std::cerr << std::endl;
+ std::cerr << "usage:\tbamtools merge [--out FILE] --in <BAM file1> [BAM file2] [BAM file3]..." << std::endl;
+ std::cerr << "\t--in\tInput BAM file(s)\t\t[at least 1 req'd]" << std::endl;
+ std::cerr << "\t--out\tDestination for merge results\t[default=stdout]" << std::endl;
+ std::cerr << std::endl;
+ return 0;
+}
+
+int RunBamMerge(int argc, char* argv[]) {
+
+ // only 'bamtool merge', show help
+ if ( argc == 2 ) return BamMergeHelp();
+
+ // else parse command line for args
+ GetOpt options(argc, argv, 1);
+
+ std::string outputFilename = "";
+ options.addOption('o', "out", &outputFilename);
+
+ std::vector<std::string> inputFilenames;
+ options.addVariableLengthOption("in", &inputFilenames);
+
+ if ( !options.parse() || inputFilenames.empty() ) return BamMergeHelp();
+ if ( outputFilename.empty() ) { outputFilename = "stdout"; }
+
+ // opens the BAM files without checking for indexes
+ BamMultiReader reader;
+ reader.Open(inputFilenames, false);
+
+ // retrieve header & reference dictionary info
+ std::string mergedHeader = reader.GetHeaderText();
+ RefVector references = reader.GetReferenceData();
+
+ // open BamWriter
+ BamWriter writer;
+ writer.Open(outputFilename, mergedHeader, references);
+
+ // store alignments to output file
+ BamAlignment bAlignment;
+ while (reader.GetNextAlignment(bAlignment)) {
+ writer.SaveAlignment(bAlignment);
+ }
+
+ // clean & exit
+ reader.Close();
+ writer.Close();
+ return 0;
+}
+
+} // namespace BamTools
+
+#endif // BAMTOOLS_MERGE_H
\ No newline at end of file
--- /dev/null
+// ***************************************************************************
+// bamtools_sam.h (c) 2010 Derek Barnett, Erik Garrison
+// Marth Lab, Department of Biology, Boston College
+// All rights reserved.
+// ---------------------------------------------------------------------------
+// Last modified: 26 May 2010
+// ---------------------------------------------------------------------------
+// Prints a BAM file in the text-based SAM format.
+// ***************************************************************************
+
+#ifndef BAMTOOLS_SAM_H
+#define BAMTOOLS_SAM_H
+
+#include <cstdlib>
+#include <iostream>
+#include <string>
+
+#include "BamReader.h"
+#include "bamtools_getopt.h"
+
+namespace BamTools {
+
+int BamSamHelp(void) {
+
+ // '--head' makes more sense than '--num' from a Unix perspective, but could be confusing with header info ??
+ // but this is also only the default case (from the beginning of the file)
+ // do we want to add a region specifier, eg 'chr2:1000..1500'? In this case, '--num' still makes sense (give me up to N alignments from this region)
+
+ std::cerr << std::endl;
+ std::cerr << "usage:\tbamtools sam [--in BAM file] [--num N] [--no_header]" << std::endl;
+ std::cerr << "\t-i, --in\tInput BAM file to generate SAM-format\t\t\t[default=stdin]" << std::endl;
+ std::cerr << "\t-n, --num N\tOnly print up to N alignments from beginning of file\t\t[default=50*]" << endl;
+ std::cerr << "\t--no_header\tOmits SAM header information from output (alignments only)\t[default=off]" << std::endl;
+ std::cerr << std::endl;
+ std::cerr << "* - By default bamtools sam will print all alignments in SAM format." << std::endl;
+ std::cerr << " However if '-n' or '--num' is included with no N, the default of 50 is used." << std::endl;
+ std::cerr << std::endl;
+ return 0;
+}
+
+static RefVector references;
+
+void PrintSAM(const BamAlignment& a) {
+
+ // tab-delimited
+ // <QNAME> <FLAG> <RNAME> <POS> <MAPQ> <CIGAR> <MRNM> <MPOS> <ISIZE> <SEQ> <QUAL> [ <TAG>:<VTYPE>:<VALUE> [...] ]
+
+ // ******************************* //
+ // ** NOT FULLY IMPLEMENTED YET ** //
+ //******************************** //
+ //
+ // Todo : build CIGAR string
+ // build TAG string
+ // there are some quirks, per the spec, regarding when to use '=' or not
+ //
+ // ******************************* //
+
+ //
+ // do validity check on RefID / MateRefID ??
+ //
+
+ // build CIGAR string
+ std::string cigarString("CIGAR:NOT YET");
+
+ // build TAG string
+ std::string tagString("TAG:NOT YET");
+
+ // print BamAlignment to stdout in SAM format
+ std::cout << a.Name << '\t'
+ << a.AlignmentFlag << '\t'
+ << references[a.RefID].RefName << '\t'
+ << a.Position << '\t'
+ << a.MapQuality << '\t'
+ << cigarString << '\t'
+ << ( a.IsPaired() ? references[a.MateRefID].RefName : "*" ) << '\t'
+ << ( a.IsPaired() ? a.MatePosition : 0 ) << '\t'
+ << ( a.IsPaired() ? a.InsertSize : 0 ) << '\t'
+ << a.QueryBases << '\t'
+ << a.Qualities << '\t'
+ << tagString << std::endl;
+}
+
+int RunBamSam(int argc, char* argv[]) {
+
+ // else parse command line for args
+ GetOpt options(argc, argv, 1);
+
+ std::string inputFilename;
+ options.addOption('i', "in", &inputFilename);
+
+ std::string numberString;
+ options.addOptionalOption('n', "num", &numberString, "50");
+
+ bool isOmittingHeader;
+ options.addSwitch("no_header", &isOmittingHeader);
+
+ if ( !options.parse() ) return BamCoverageHelp();
+ if ( inputFilename.empty() ) { inputFilename = "stdin"; }
+
+ // maxNumberOfAlignments = all (if nothing specified)
+ // = 50 (if '-n' or '--num' but no N)
+ // = N (if '-n N' or '--num N')
+ int maxNumberOfAlignments = -1;
+ if ( !numberString.empty() ) { maxNumberOfAlignments = atoi(numberString.c_str()); }
+
+ // open our BAM reader
+ BamReader reader;
+ reader.Open(inputFilename);
+
+ // if header desired, retrieve and print to stdout
+ if ( !isOmittingHeader ) {
+ std::string header = reader.GetHeaderText();
+ std::cout << header << std::endl;
+ }
+
+ // store reference data
+ references = reader.GetReferenceData();
+
+ // print all alignments to stdout in SAM format
+ if ( maxNumberOfAlignments < 0 ) {
+ BamAlignment ba;
+ while( reader.GetNextAlignment(ba) ) {
+ PrintSAM(ba);
+ }
+ }
+
+ // print first N alignments to stdout in SAM format
+ else {
+ BamAlignment ba;
+ int alignmentsPrinted = 0;
+ while ( reader.GetNextAlignment(ba) && (alignmentsPrinted < maxNumberOfAlignments) ) {
+ PrintSAM(ba);
+ ++alignmentsPrinted;
+ }
+ }
+
+ // clean & exit
+ reader.Close();
+ return 0;
+}
+
+} // namespace BamTools
+
+#endif // BAMTOOLS_SAM_H
\ No newline at end of file
--- /dev/null
+// ***************************************************************************
+// bamtools_sortt.h (c) 2010 Derek Barnett, Erik Garrison
+// Marth Lab, Department of Biology, Boston College
+// All rights reserved.
+// ---------------------------------------------------------------------------
+// Last modified: 26 May 2010
+// ---------------------------------------------------------------------------
+// Sorts an input BAM file (default by position) and stores in a new BAM file.
+// ***************************************************************************
+
+#ifndef BAMTOOLS_SORT_H
+#define BAMTOOLS_SORT_H
+
+#include <iostream>
+#include <string>
+
+#include "BamReader.h"
+#include "bamtools_getopt.h"
+
+namespace BamTools {
+
+int BamSortHelp(void) {
+ std::cerr << std::endl;
+ std::cerr << "usage:\tbamtools sort [--in BAM file] [--out sorted BAM file]" << std::endl;
+ std::cerr << "\t-i, --in\tInput BAM file to sort\t[default=stdin]" << std::endl;
+ std::cerr << "\t-o. --out\tDestination of sorted BAM file\t[default=stdout]" << std::endl;
+ std::cerr << std::endl;
+ return 0;
+}
+
+int RunBamSort(int argc, char* argv[]) {
+
+ // else parse command line for args
+ GetOpt options(argc, argv, 1);
+
+ std::string inputFilename;
+ options.addOption('i', "in", &inputFilename);
+
+ std::string outputFilename;
+ options.addOption('o', "out", &outputFilename);
+
+ if ( !options.parse() ) return BamCoverageHelp();
+ if ( inputFilename.empty() ) { inputFilename = "stdin"; }
+ if ( outputFilename.empty() ) { outputFilename = "stdout"; }
+
+ // open our BAM reader
+// BamReader reader;
+// reader.Open(inputFilename);
+//
+// // retrieve header & reference dictionary info
+// std::string header = reader.GetHeaderText();
+// RefVector references = reader.GetReferenceData();
+//
+// BamWriter writer;
+// writer.Open(outputFilename, header, references);
+//
+ // sort BAM file
+ std::cerr << "Sorting " << inputFilename << std::endl;
+ std::cerr << "Saving sorted BAM in " << outputFilename << endl;
+ std::cerr << "FEATURE NOT YET IMPLEMENTED!" << std::endl;
+
+ // clean & exit
+// reader.Close();
+// writer.Close();
+ return 0;
+}
+
+} // namespace BamTools
+
+#endif // BAMTOOLS_SORT_H
\ No newline at end of file
--- /dev/null
+// ***************************************************************************
+// bamtools_stats.h (c) 2010 Derek Barnett, Erik Garrison
+// Marth Lab, Department of Biology, Boston College
+// All rights reserved.
+// ---------------------------------------------------------------------------
+// Last modified: 26 May 2010
+// ---------------------------------------------------------------------------
+// Prints general statistics for a single BAM file
+//
+// ** Expand to multiple??
+//
+// ***************************************************************************
+
+#ifndef BAMTOOLS_STATS_H
+#define BAMTOOLS_STATS_H
+
+#include <iostream>
+#include <string>
+
+#include "BamReader.h"
+#include "bamtools_getopt.h"
+
+namespace BamTools {
+
+int BamStatsHelp(void) {
+ std::cerr << std::endl;
+ std::cerr << "usage:\tbamtools stats [--in BAM file]" << std::endl;
+ std::cerr << "\t-i, --in\tInput BAM file to calculate general stats\t[default=stdin]" << std::endl;
+ std::cerr << std::endl;
+ return 0;
+}
+
+int RunBamStats(int argc, char* argv[]) {
+
+ // else parse command line for args
+ GetOpt options(argc, argv, 1);
+
+ std::string inputFilename;
+ options.addOption('t', "inp", &inputFilename);
+
+ if ( !options.parse() ) return BamStatsHelp();
+ if ( inputFilename.empty() ) { inputFilename = "stdin"; }
+
+ // open our BAM reader
+// BamReader reader;
+// reader.Open(inputFilename);
+
+ // calculate general stats
+ std::cerr << "Calculating general stats for " << inputFilename << std::endl;
+ std::cerr << "FEATURE NOT YET IMPLEMENTED!" << std::endl;
+
+ // clean & exit
+// reader.Close();
+ return 0;
+}
+
+} // namespace BamTools
+
+#endif // BAMTOOLS_STATS_H
\ No newline at end of file