#include "bamtools_options.h"
#include "bamtools_utilities.h"
#include "BamReader.h"
+#include "BamMultiReader.h"
using namespace std;
using namespace BamTools;
struct CountTool::CountSettings {
// flags
- bool HasBamIndexFilename;
- bool HasInputBamFilename;
+ bool HasInput;
bool HasRegion;
// filenames
- string BamIndexFilename;
- string InputBamFilename;
+ vector<string> InputFiles;
string Region;
// constructor
CountSettings(void)
- : HasBamIndexFilename(false)
- , HasInputBamFilename(false)
+ : HasInput(false)
, HasRegion(false)
- , InputBamFilename(Options::StandardIn())
{ }
};
, m_settings(new CountSettings)
{
// set program details
- Options::SetProgramInfo("bamtools count", "prints alignment counts for a BAM file", "-in <filename> [-region <REGION> [-index <filename>]]");
+ Options::SetProgramInfo("bamtools count", "prints alignment counts for a BAM file", "-in <filename> [-region <REGION>]");
// set up options
OptionGroup* IO_Opts = Options::CreateOptionGroup("Input & Output");
- Options::AddValueOption("-in", "BAM filename", "the input BAM file", "", m_settings->HasInputBamFilename, m_settings->InputBamFilename, IO_Opts, Options::StandardIn());
- Options::AddValueOption("-index", "BAM index filename", "the BAM index file", "", m_settings->HasBamIndexFilename, m_settings->BamIndexFilename, IO_Opts);
+ Options::AddValueOption("-in", "BAM filename", "the input BAM file(s)", "", m_settings->HasInput, m_settings->InputFiles, IO_Opts);
+ //Options::AddValueOption("-index", "BAM index filename", "the BAM index file", "", m_settings->HasBamIndexFilename, m_settings->BamIndexFilename, IO_Opts);
OptionGroup* FilterOpts = Options::CreateOptionGroup("Filters");
- Options::AddValueOption("-region", "REGION", "genomic region. Index file is recommended for better performance. See \'bamtools help index\' for more details on creating one", "", m_settings->HasRegion, m_settings->Region, FilterOpts);
+ Options::AddValueOption("-region", "REGION", "genomic region. Index file is recommended for better performance, and is read automatically if it exists as <filename>.bai or <filename>.bti. See \'bamtools help index\' for more details on creating one", "", m_settings->HasRegion, m_settings->Region, FilterOpts);
}
CountTool::~CountTool(void) {
// parse command line arguments
Options::Parse(argc, argv, 1);
- //open our BAM reader
- BamReader reader;
- reader.Open(m_settings->InputBamFilename);
+ if ( !m_settings->HasInput )
+ m_settings->InputFiles.push_back(Options::StandardIn());
+
+ BamMultiReader reader;
+ reader.Open(m_settings->InputFiles, false, true);
// alignment counter
int alignmentCount(0);
// if no region specified, count entire file
if ( !m_settings->HasRegion ) {
BamAlignment al;
- while ( reader.GetNextAlignment(al) )
+ while ( reader.GetNextAlignmentCore(al) )
++alignmentCount;
- }
+ }
// more complicated - region specified
else {
- Region region;
+ BamRegion region;
if ( Utilities::ParseRegionString(m_settings->Region, reader, region) ) {
- // has index, we can jump directly to
- if ( m_settings->HasBamIndexFilename ) {
-
- // this is kind of a hack...?
- // re-open BamReader, this time with the index file
- reader.Close();
- reader.Open(m_settings->InputBamFilename, m_settings->BamIndexFilename);
+ // check if there are index files *.bai/*.bti corresponding to the input files
+ bool hasDefaultIndex = false;
+ bool hasBamtoolsIndex = false;
+ bool hasNoIndex = false;
+ int defaultIndexCount = 0;
+ int bamtoolsIndexCount = 0;
+ for (vector<string>::const_iterator f = m_settings->InputFiles.begin(); f != m_settings->InputFiles.end(); ++f) {
- // attempt Jump(), catch error
- // if ( !reader.Jump(region.StartChromID, region.StartPosition) ) {
- // foundError = true;
- // errorStream << "Could not jump to desired REGION: " << m_settings->Region << endl;
- // }
- //
- if ( !reader.SetRegion(region.StartChromID, region.StartPosition, region.StopChromID, region.StopPosition) ) {
- foundError = true;
- errorStream << "Could not set BamReader region to REGION: " << m_settings->Region << endl;
+ if ( Utilities::FileExists(*f + ".bai") ) {
+ hasDefaultIndex = true;
+ ++defaultIndexCount;
+ }
+
+ if ( Utilities::FileExists(*f + ".bti") ) {
+ hasBamtoolsIndex = true;
+ ++bamtoolsIndexCount;
+ }
+
+ if ( !hasDefaultIndex && !hasBamtoolsIndex ) {
+ hasNoIndex = true;
+ cerr << "*WARNING - could not find index file for " << *f
+ << ", parsing whole file(s) to get alignment counts for target region"
+ << " (could be slow)" << endl;
+ break;
}
-
}
- else {
-
- // read through sequentially, until first overlapping read is found
+ // determine if index file types are heterogeneous
+ bool hasDifferentIndexTypes = false;
+ if ( defaultIndexCount > 0 && bamtoolsIndexCount > 0 ) {
+ hasDifferentIndexTypes = true;
+ cerr << "*WARNING - different index file formats found"
+ << ", parsing whole file(s) to get alignment counts for target region"
+ << " (could be slow)" << endl;
+ }
+
+ // if any input file has no index, or if input files use different index formats
+ // can't use BamMultiReader to jump directly (**for now**)
+ if ( hasNoIndex || hasDifferentIndexTypes ) {
+
+ // read through sequentially, counting all overlapping reads
BamAlignment al;
- bool alignmentFound(false);
- while( reader.GetNextAlignment(al) ) {
- if ( (al.RefID == region.StartChromID ) && ( (al.Position + al.Length) >= region.StartPosition) ) {
- alignmentFound = true;
- break;
+ while( reader.GetNextAlignmentCore(al) ) {
+ if ( (al.RefID >= region.LeftRefID) && ( (al.Position + al.Length) >= region.LeftPosition ) &&
+ (al.RefID <= region.RightRefID) && ( al.Position <= region.RightPosition) )
+ {
+ ++alignmentCount;
}
}
-
- // if overlapping alignment found (not EOF), increment counter
- if ( alignmentFound) ++ alignmentCount;
}
- // -----------------------------
- // count alignments until stop hit
-
- if ( !foundError ) {
- // while valid alignment AND
- // ( either current ref is before stopRefID OR
- // current ref stopRefID but we're before stopPos )
- BamAlignment al;
- //while ( reader.GetNextAlignment(al) && ((al.RefID < region.StopChromID ) || (al.RefID == region.StopChromID && al.Position <= region.StopPosition)) )
- // ++alignmentCount;
- //
- while ( reader.GetNextAlignment(al) )
- ++alignmentCount;
-
+ // has index file for each input file (and same format)
+ else {
+
+ // this is kind of a hack...?
+ BamMultiReader reader;
+ reader.Open(m_settings->InputFiles, true, true, hasDefaultIndex );
+
+ if ( !reader.SetRegion(region.LeftRefID, region.LeftPosition, region.RightRefID, region.RightPosition) ) {
+ foundError = true;
+ errorStream << "Could not set BamReader region to REGION: " << m_settings->Region << endl;
+ } else {
+ BamAlignment al;
+ while ( reader.GetNextAlignmentCore(al) )
+ ++alignmentCount;
+ }
}
-
+
} else {
foundError = true;
errorStream << "Could not parse REGION: " << m_settings->Region << endl;
projects / bamtools.git / blobdiff