1 // ***************************************************************************
2 // bamtools_coverage.cpp (c) 2010 Derek Barnett, Erik Garrison
3 // Marth Lab, Department of Biology, Boston College
4 // ---------------------------------------------------------------------------
5 // Last modified: 7 April 2011
6 // ---------------------------------------------------------------------------
7 // Prints coverage data for a single BAM file
8 // ***************************************************************************
10 #include "bamtools_coverage.h"
12 #include <api/BamReader.h>
13 #include <utils/bamtools_options.h>
14 #include <utils/bamtools_pileup_engine.h>
15 using namespace BamTools;
25 // ---------------------------------------------
26 // CoverageVisitor implementation
28 class CoverageVisitor : public PileupVisitor {
31 CoverageVisitor(const RefVector& references, ostream* out)
33 , m_references(references)
36 ~CoverageVisitor(void) { }
38 // PileupVisitor interface implementation
40 // prints coverage results ( tab-delimited )
41 void Visit(const PileupPosition& pileupData) {
42 *m_out << m_references[pileupData.RefId].RefName << "\t"
43 << pileupData.Position << "\t"
44 << pileupData.PileupAlignments.size() << endl;
48 RefVector m_references;
52 } // namespace BamTools
54 // ---------------------------------------------
55 // CoverageSettings implementation
57 struct CoverageTool::CoverageSettings {
64 string InputBamFilename;
65 string OutputFilename;
68 CoverageSettings(void)
70 , HasOutputFile(false)
71 , InputBamFilename(Options::StandardIn())
72 , OutputFilename(Options::StandardOut())
76 // ---------------------------------------------
77 // CoverageToolPrivate implementation
79 struct CoverageTool::CoverageToolPrivate {
83 CoverageToolPrivate(CoverageTool::CoverageSettings* settings)
84 : m_settings(settings)
88 ~CoverageToolPrivate(void) { }
96 CoverageTool::CoverageSettings* m_settings;
98 RefVector m_references;
101 bool CoverageTool::CoverageToolPrivate::Run(void) {
103 // if output filename given
105 if ( m_settings->HasOutputFile ) {
107 // open output file stream
108 outFile.open(m_settings->OutputFilename.c_str());
110 cerr << "bamtools coverage ERROR: could not open " << m_settings->OutputFilename
111 << " for output" << endl;
115 // set m_out to file's streambuf
116 m_out.rdbuf(outFile.rdbuf());
119 //open our BAM reader
121 if ( !reader.Open(m_settings->InputBamFilename) ) {
122 cerr << "bamtools coverage ERROR: could not open input BAM file: " << m_settings->InputBamFilename << endl;
126 // retrieve references
127 m_references = reader.GetReferenceData();
129 // set up our output 'visitor'
130 CoverageVisitor* cv = new CoverageVisitor(m_references, &m_out);
132 // set up pileup engine with 'visitor'
134 pileup.AddVisitor(cv);
136 // process input data
138 while ( reader.GetNextAlignment(al) )
139 pileup.AddAlignment(al);
143 if ( m_settings->HasOutputFile )
152 // ---------------------------------------------
153 // CoverageTool implementation
155 CoverageTool::CoverageTool(void)
157 , m_settings(new CoverageSettings)
160 // set program details
161 Options::SetProgramInfo("bamtools coverage", "prints coverage data for a single BAM file", "[-in <filename>] [-out <filename>]");
164 OptionGroup* IO_Opts = Options::CreateOptionGroup("Input & Output");
165 Options::AddValueOption("-in", "BAM filename", "the input BAM file", "", m_settings->HasInputFile, m_settings->InputBamFilename, IO_Opts, Options::StandardIn());
166 Options::AddValueOption("-out", "filename", "the output file", "", m_settings->HasOutputFile, m_settings->OutputFilename, IO_Opts, Options::StandardOut());
169 CoverageTool::~CoverageTool(void) {
178 int CoverageTool::Help(void) {
179 Options::DisplayHelp();
183 int CoverageTool::Run(int argc, char* argv[]) {
185 // parse command line arguments
186 Options::Parse(argc, argv, 1);
188 // initialize CoverageTool with settings
189 m_impl = new CoverageToolPrivate(m_settings);
191 // run CoverageTool, return success/fail