From 889210eae78b7e1a186b67aa0846c038f03906ab Mon Sep 17 00:00:00 2001 From: derek Date: Tue, 10 May 2011 19:28:31 -0400 Subject: [PATCH] Fixed minimum offset error in BAI jumping --- src/api/internal/BamReader_p.cpp | 16 ++++++++++++---- src/api/internal/BamStandardIndex_p.cpp | 15 +++++++++++---- src/api/internal/BamWriter_p.cpp | 4 ++-- 3 files changed, 25 insertions(+), 10 deletions(-) diff --git a/src/api/internal/BamReader_p.cpp b/src/api/internal/BamReader_p.cpp index c2afb4c..5daa1bf 100644 --- a/src/api/internal/BamReader_p.cpp +++ b/src/api/internal/BamReader_p.cpp @@ -3,7 +3,7 @@ // Marth Lab, Department of Biology, Boston College // All rights reserved. // --------------------------------------------------------------------------- -// Last modified: 5 April 2011 (DB) +// Last modified: 10 May 2011 (DB) // --------------------------------------------------------------------------- // Provides the basic functionality for reading BAM files // *************************************************************************** @@ -306,16 +306,24 @@ bool BamReaderPrivate::Open(const string& filename) { Close(); // attempt to open BgzfStream for reading - if ( !m_stream.Open(filename, "rb") ) + if ( !m_stream.Open(filename, "rb") ) { + cerr << "BamReader ERROR: Could not open BGZF stream for " << filename << endl; return false; + } // attempt to load header data - if ( !LoadHeaderData() ) + if ( !LoadHeaderData() ) { + cerr << "BamReader ERROR: Could not load header data for " << filename << endl; + Close(); return false; + } // attempt to load reference data - if ( !LoadReferenceData() ) + if ( !LoadReferenceData() ) { + cerr << "BamReader ERROR: Could not load reference data for " << filename << endl; + Close(); return false; + } // if all OK, store filename & offset of first alignment m_filename = filename; diff --git a/src/api/internal/BamStandardIndex_p.cpp b/src/api/internal/BamStandardIndex_p.cpp index d29c3c5..c148d5d 100644 --- a/src/api/internal/BamStandardIndex_p.cpp +++ b/src/api/internal/BamStandardIndex_p.cpp @@ -3,7 +3,7 @@ // Marth Lab, Department of Biology, Boston College // All rights reserved. // --------------------------------------------------------------------------- -// Last modified: 27 April 2011 (DB) +// Last modified: 10 May 2011 (DB) // --------------------------------------------------------------------------- // Provides index operations for the standardized BAM index format (".bai") // *************************************************************************** @@ -136,7 +136,7 @@ bool BamStandardIndex::CalculateCandidateOffsets(const BaiReferenceSummary& refS // store alignment chunk's start offset // if its stop offset is larger than our 'minOffset' - if ( chunkStop > minOffset ) + if ( chunkStop >= minOffset ) offsets.push_back(chunkStart); } @@ -414,8 +414,10 @@ bool BamStandardIndex::GetOffsets(const BamRegion& region, vector& offs // set up region boundaries based on actual BamReader data uint32_t begin; uint32_t end; - if ( !AdjustRegion(region, begin, end) ) + if ( !AdjustRegion(region, begin, end) ) { + cerr << "BamStandardIndex ERROR: cannot calculate offsets on invalid region" << endl; return false; + } // retrieve all candidate bin IDs for region set candidateBins; @@ -426,8 +428,10 @@ bool BamStandardIndex::GetOffsets(const BamRegion& region, vector& offs const uint64_t& minOffset = CalculateMinOffset(refSummary, begin); // attempt to use reference summary, minOffset, & candidateBins to calculate offsets - if ( !CalculateCandidateOffsets(refSummary, minOffset, candidateBins, offsets) ) + if ( !CalculateCandidateOffsets(refSummary, minOffset, candidateBins, offsets) ) { + cerr << "Could not caluclate candidate offsets for region" << endl; return false; + } // ensure that offsets are sorted before returning sort( offsets.begin(), offsets.end() ); @@ -684,6 +688,8 @@ void BamStandardIndex::SaveAlignmentChunkToBin(BaiBinMap& binMap, // create new alignment chunk BaiAlignmentChunk newChunk(currentOffset, lastOffset); + + // if no entry exists yet for this bin, create one and store alignment chunk BaiBinMap::iterator binIter = binMap.find(currentBin); if ( binIter == binMap.end() ) { @@ -822,6 +828,7 @@ bool BamStandardIndex::SummarizeLinearOffsets(BaiReferenceSummary& refSummary) { } bool BamStandardIndex::SummarizeReference(BaiReferenceSummary& refSummary) { + bool loadedOk = true; loadedOk &= SummarizeBins(refSummary); loadedOk &= SummarizeLinearOffsets(refSummary); diff --git a/src/api/internal/BamWriter_p.cpp b/src/api/internal/BamWriter_p.cpp index 8f6c30d..fbe64db 100644 --- a/src/api/internal/BamWriter_p.cpp +++ b/src/api/internal/BamWriter_p.cpp @@ -3,7 +3,7 @@ // Marth Lab, Department of Biology, Boston College // All rights reserved. // --------------------------------------------------------------------------- -// Last modified: 19 April 2011 (DB) +// Last modified: 10 May 2011 (DB) // --------------------------------------------------------------------------- // Provides the basic functionality for producing BAM files // *************************************************************************** @@ -258,7 +258,7 @@ void BamWriterPrivate::SaveAlignment(const BamAlignment& al) { // write the base qualities char* pBaseQualities = (char*)al.Qualities.data(); - for ( unsigned int i = 0; i < queryLength; i++ ) + for ( unsigned int i = 0; i < queryLength; ++i ) pBaseQualities[i] -= 33; // FASTQ conversion m_stream.Write(pBaseQualities, queryLength); -- 2.39.2