// Marth Lab, Department of Biology, Boston College\r
// All rights reserved.\r
// ---------------------------------------------------------------------------\r
-// Last modified: 7 September 2010 (DB)\r
+// Last modified: 10 September 2010 (DB)\r
// ---------------------------------------------------------------------------\r
// Uses BGZF routines were adapted from the bgzf.c code developed at the Broad\r
// Institute.\r
// general file data\r
BgzfData mBGZF;\r
string HeaderText;\r
- //BamIndex Index;\r
- BamIndex* NewIndex;\r
+ BamIndex* Index;\r
RefVector References;\r
bool IsIndexLoaded;\r
int64_t AlignmentsBeginOffset;\r
bool IsBigEndian;\r
\r
// user-specified region values\r
- BamRegion Region;\r
- bool IsLeftBoundSpecified;\r
- bool IsRightBoundSpecified;\r
- \r
- bool IsRegionSpecified;\r
+ BamRegion Region; \r
int CurrentRefID;\r
int CurrentLeft;\r
\r
\r
// file operations\r
void Close(void);\r
- bool Jump(int refID, int position);\r
bool Open(const std::string& filename, \r
const std::string& indexFilename, \r
const bool lookForIndex, \r
void BamReader::Close(void) { d->Close(); }\r
bool BamReader::IsIndexLoaded(void) const { return d->IsIndexLoaded; }\r
bool BamReader::IsOpen(void) const { return d->mBGZF.IsOpen; }\r
-bool BamReader::Jump(int refID, int position) \r
-{ \r
- d->Region.LeftRefID = refID;\r
- d->Region.LeftPosition = position;\r
- d->IsLeftBoundSpecified = true;\r
- d->IsRightBoundSpecified = false;\r
- return d->Jump(refID, position); \r
-}\r
+bool BamReader::Jump(int refID, int position) { return d->SetRegion( BamRegion(refID, position) ); }\r
bool BamReader::Open(const std::string& filename, \r
const std::string& indexFilename, \r
const bool lookForIndex, \r
\r
// constructor\r
BamReader::BamReaderPrivate::BamReaderPrivate(BamReader* parent)\r
- : NewIndex(0)\r
+ : Index(0)\r
, IsIndexLoaded(false)\r
, AlignmentsBeginOffset(0)\r
- , IsLeftBoundSpecified(false)\r
- , IsRightBoundSpecified(false)\r
- , IsRegionSpecified(false)\r
, CurrentRefID(0)\r
, CurrentLeft(0)\r
, Parent(parent)\r
\r
// clear index data structure\r
void BamReader::BamReaderPrivate::ClearIndex(void) {\r
- delete NewIndex;\r
- NewIndex = 0;\r
+ delete Index;\r
+ Index = 0;\r
IsIndexLoaded = false;\r
}\r
\r
HeaderText.clear();\r
\r
// clear out region flags\r
- IsLeftBoundSpecified = false;\r
- IsRightBoundSpecified = false;\r
- IsRegionSpecified = false;\r
+ Region.clear();\r
}\r
\r
// creates index for BAM file, saves to file\r
\r
// create index based on type requested\r
if ( useStandardIndex ) \r
- NewIndex = new BamStandardIndex(&mBGZF, Parent, IsBigEndian);\r
+ Index = new BamStandardIndex(&mBGZF, Parent, IsBigEndian);\r
// create BamTools 'custom' index\r
else\r
- NewIndex = new BamToolsIndex(&mBGZF, Parent, IsBigEndian);\r
+ Index = new BamToolsIndex(&mBGZF, Parent, IsBigEndian);\r
\r
// build new index\r
bool ok = true;\r
- ok &= NewIndex->Build();\r
+ ok &= Index->Build();\r
IsIndexLoaded = ok;\r
\r
// attempt to save index data to file\r
- ok &= NewIndex->Write(Filename); \r
+ ok &= Index->Write(Filename); \r
\r
// return success/fail of both building & writing index\r
return ok;\r
bAlignment.SupportData.HasCoreOnly = true;\r
\r
// if region not specified, return success\r
- if ( !IsLeftBoundSpecified ) return true;\r
+ if ( !Region.isLeftBoundSpecified() ) return true;\r
\r
// determine region state (before, within, after)\r
BamReader::BamReaderPrivate::RegionState state = IsOverlap(bAlignment);\r
// check alignment start against right bound cutoff\r
\r
// if full region of interest was given\r
- if ( IsRightBoundSpecified ) {\r
+ if ( Region.isRightBoundSpecified() ) {\r
\r
// read starts on right bound reference, but AFTER right bound position\r
if ( bAlignment.RefID == Region.RightRefID && bAlignment.Position > Region.RightPosition )\r
return BEFORE_REGION;\r
}\r
\r
-// jumps to specified region(refID, leftBound) in BAM file, returns success/fail\r
-bool BamReader::BamReaderPrivate::Jump(int refID, int position) {\r
-\r
- // -----------------------------------------------------------------------\r
- // check for existing index \r
- if ( !IsIndexLoaded || NewIndex == 0 ) return false; \r
- // see if reference has alignments\r
- if ( !NewIndex->HasAlignments(refID) ) return false; \r
- // make sure position is valid\r
- if ( position > References.at(refID).RefLength ) return false;\r
- \r
- // determine possible offsets\r
- vector<int64_t> offsets;\r
- if ( !NewIndex->GetOffsets(Region, IsRightBoundSpecified, offsets) ) {\r
- fprintf(stderr, "ERROR: Could not jump: unable to calculate offset for specified region.\n");\r
- return false;\r
- }\r
- \r
- // iterate through offsets\r
- BamAlignment bAlignment;\r
- bool result = true;\r
- for ( vector<int64_t>::const_iterator o = offsets.begin(); o != offsets.end(); ++o) {\r
- \r
- // attempt seek & load first available alignment\r
- result &= mBGZF.Seek(*o);\r
- LoadNextAlignment(bAlignment);\r
- \r
- // if this alignment corresponds to desired position\r
- // return success of seeking back to 'current offset'\r
- if ( (bAlignment.RefID == refID && bAlignment.Position + bAlignment.Length > position) || (bAlignment.RefID > refID) ) {\r
- if ( o != offsets.begin() ) --o;\r
- return mBGZF.Seek(*o);\r
- }\r
- }\r
- \r
- return result;\r
-}\r
-\r
// load BAM header data\r
void BamReader::BamReaderPrivate::LoadHeaderData(void) {\r
\r
\r
// attempt to load BamIndex based on current Filename provided & preferStandardIndex flag\r
const BamIndex::PreferredIndexType type = (preferStandardIndex ? BamIndex::STANDARD : BamIndex::BAMTOOLS);\r
- NewIndex = BamIndex::FromBamFilename(Filename, &mBGZF, Parent, IsBigEndian, type);\r
+ Index = BamIndex::FromBamFilename(Filename, &mBGZF, Parent, IsBigEndian, type);\r
\r
// if null, return failure\r
- if ( NewIndex == 0 ) return false;\r
+ if ( Index == 0 ) return false;\r
\r
// generate proper IndexFilename based on type of index created\r
- IndexFilename = Filename + NewIndex->Extension();\r
+ IndexFilename = Filename + Index->Extension();\r
}\r
\r
else {\r
// attempt to load BamIndex based on IndexFilename provided by client\r
- NewIndex = BamIndex::FromIndexFilename(IndexFilename, &mBGZF, Parent, IsBigEndian);\r
+ Index = BamIndex::FromIndexFilename(IndexFilename, &mBGZF, Parent, IsBigEndian);\r
\r
// if null, return failure\r
- if ( NewIndex == 0 ) return false; \r
+ if ( Index == 0 ) return false; \r
}\r
\r
// an index file was found\r
// return success of loading the index data from file\r
- IsIndexLoaded = NewIndex->Load(IndexFilename);\r
+ IsIndexLoaded = Index->Load(IndexFilename);\r
return IsIndexLoaded;\r
}\r
\r
// returns BAM file pointer to beginning of alignment data\r
bool BamReader::BamReaderPrivate::Rewind(void) {\r
\r
- // rewind to first alignment\r
+ // rewind to first alignment, return false if unable to seek\r
if ( !mBGZF.Seek(AlignmentsBeginOffset) ) return false;\r
\r
- // retrieve first alignment data\r
+ // retrieve first alignment data, return false if unable to read\r
BamAlignment al;\r
if ( !LoadNextAlignment(al) ) return false;\r
\r
// reset default region info using first alignment in file\r
- Region.LeftRefID = al.RefID;\r
- Region.LeftPosition = al.Position;\r
- Region.RightRefID = -1;\r
- Region.RightPosition = -1;\r
- IsLeftBoundSpecified = false;\r
- IsRightBoundSpecified = false; \r
-\r
- // rewind back to before first alignment\r
+ Region.clear();\r
+ Region.LeftRefID = al.RefID;\r
+ Region.LeftPosition = al.Position;\r
+\r
+ // rewind back to beginning of first alignment\r
// return success/fail of seek\r
return mBGZF.Seek(AlignmentsBeginOffset);\r
}\r
\r
-// sets a region of interest (with left & right bound reference/position)\r
-// attempts a Jump() to left bound as well\r
-// returns success/failure of Jump()\r
+// asks Index to attempt a Jump() to specified region\r
+// returns success/failure\r
bool BamReader::BamReaderPrivate::SetRegion(const BamRegion& region) {\r
+ \r
+ // clear out any prior BamReader region data\r
+ //\r
+ // N.B. - this is cleared so that BamIndex now has free reign to call\r
+ // GetNextAlignmentCore() and do overlap checking without worrying about BamReader \r
+ // performing any overlap checking of its own and moving on to the next read... Calls \r
+ // to GetNextAlignmentCore() with no Region set, simply return the next alignment.\r
+ // This ensures that the Index is able to do just that. (All without exposing \r
+ // LoadNextAlignment() to the public API, and potentially confusing clients with the nomenclature)\r
+ Region.clear();\r
+ \r
+ // check for existing index \r
+ if ( !IsIndexLoaded || Index == 0 ) return false; \r
\r
- // save region of interest\r
+ // attempt jump to user-specified region, return false if failed\r
+ if ( !Index->Jump(region) ) return false;\r
+ \r
+ // if successful, save region data locally, return success\r
Region = region;\r
- \r
- // set flags\r
- if ( region.LeftRefID >= 0 && region.LeftPosition >= 0 ) IsLeftBoundSpecified = true;\r
- if ( region.RightRefID >= 0 && region.RightPosition >= 0 ) IsRightBoundSpecified = true;\r
- \r
- // attempt jump to beginning of region, return success/fail of Jump()\r
- return Jump( Region.LeftRefID, Region.LeftPosition );\r
+ return true;\r
}\r
projects / bamtools.git / blobdiff