1 // ***************************************************************************
2 // BamRandomAccessController_p.cpp (c) 2011 Derek Barnett
3 // Marth Lab, Department of Biology, Boston College
4 // ---------------------------------------------------------------------------
5 // Last modified: 10 October 2011(DB)
6 // ---------------------------------------------------------------------------
7 // Manages random access operations in a BAM file
8 // **************************************************************************
10 #include "api/BamIndex.h"
11 #include "api/internal/BamException_p.h"
12 #include "api/internal/BamRandomAccessController_p.h"
13 #include "api/internal/BamReader_p.h"
14 #include "api/internal/BamIndexFactory_p.h"
15 using namespace BamTools;
16 using namespace BamTools::Internal;
22 BamRandomAccessController::BamRandomAccessController(void)
24 , m_hasAlignmentsInRegion(true)
27 BamRandomAccessController::~BamRandomAccessController(void) {
31 void BamRandomAccessController::AdjustRegion(const int& referenceCount) {
33 // skip if no index available
37 // see if any references in region have alignments
38 m_hasAlignmentsInRegion = false;
39 int currentId = m_region.LeftRefID;
40 const int rightBoundRefId = ( m_region.isRightBoundSpecified() ? m_region.RightRefID : referenceCount - 1 );
41 while ( currentId <= rightBoundRefId ) {
42 m_hasAlignmentsInRegion = m_index->HasAlignments(currentId);
43 if ( m_hasAlignmentsInRegion ) break;
47 // if no data found on any reference in region
48 if ( !m_hasAlignmentsInRegion )
51 // if left bound of desired region had no data, use first reference that had data
52 // otherwise, leave requested region as-is
53 if ( currentId != m_region.LeftRefID ) {
54 m_region.LeftRefID = currentId;
55 m_region.LeftPosition = 0;
59 // returns alignments' "RegionState": { Before|Overlaps|After } current region
60 BamRandomAccessController::RegionState
61 BamRandomAccessController::AlignmentState(const BamAlignment& alignment) const {
63 // if region has no left bound at all
64 if ( !m_region.isLeftBoundSpecified() )
65 return OverlapsRegion;
67 // handle unmapped reads - return AFTER region to halt processing
68 if ( alignment.RefID == -1 )
71 // if alignment is on any reference before left bound reference
72 if ( alignment.RefID < m_region.LeftRefID )
75 // if alignment is on left bound reference
76 else if ( alignment.RefID == m_region.LeftRefID ) {
78 // if alignment starts at or after left bound position
79 if ( alignment.Position >= m_region.LeftPosition) {
81 if ( m_region.isRightBoundSpecified() && // right bound is specified AND
82 m_region.LeftRefID == m_region.RightRefID && // left & right bounds on same reference AND
83 alignment.Position >= m_region.RightPosition ) // alignment starts on or after right bound position
86 // otherwise, alignment overlaps region
87 else return OverlapsRegion;
90 // alignment starts before left bound position
93 // if alignment overlaps left bound position
94 if ( alignment.GetEndPosition() > m_region.LeftPosition )
95 return OverlapsRegion;
101 // otherwise alignment is on a reference after left bound reference
104 // if region has a right bound
105 if ( m_region.isRightBoundSpecified() ) {
107 // alignment is on any reference between boundaries
108 if ( alignment.RefID < m_region.RightRefID )
109 return OverlapsRegion;
111 // alignment is on any reference after right boundary
112 else if ( alignment.RefID > m_region.RightRefID )
115 // alignment is on right bound reference
118 // if alignment starts before right bound position
119 if ( alignment.Position < m_region.RightPosition )
120 return OverlapsRegion;
126 // otherwise, alignment starts after left bound and there is no right bound given
127 else return OverlapsRegion;
131 void BamRandomAccessController::Close(void) {
136 void BamRandomAccessController::ClearIndex(void) {
143 void BamRandomAccessController::ClearRegion(void) {
145 m_hasAlignmentsInRegion = true;
148 bool BamRandomAccessController::CreateIndex(BamReaderPrivate* reader,
149 const BamIndex::IndexType& type)
151 // skip if reader is invalid
153 if ( !reader->IsOpen() ) {
154 SetErrorString("BamRandomAccessController::CreateIndex",
155 "cannot create index for unopened reader");
159 // create new index of requested type
160 BamIndex* newIndex = BamIndexFactory::CreateIndexOfType(type, reader);
161 if ( newIndex == 0 ) {
163 s << "could not create index of type: " << type;
164 SetErrorString("BamRandomAccessController::CreateIndex", s.str());
168 // attempt to build index from current BamReader file
169 if ( !newIndex->Create() ) {
170 const string indexError = newIndex->GetErrorString();
171 const string message = "could not create index: \n\t" + indexError;
172 SetErrorString("BamRandomAccessController::CreateIndex", message);
176 // save new index & return success
181 string BamRandomAccessController::GetErrorString(void) const {
182 return m_errorString;
185 bool BamRandomAccessController::HasIndex(void) const {
186 return ( m_index != 0 );
189 bool BamRandomAccessController::HasRegion(void) const {
190 return ( !m_region.isNull() );
193 bool BamRandomAccessController::IndexHasAlignmentsForReference(const int& refId) {
194 return m_index->HasAlignments(refId);
197 bool BamRandomAccessController::LocateIndex(BamReaderPrivate* reader,
198 const BamIndex::IndexType& preferredType)
200 // look up index filename, deferring to preferredType if possible
202 const string& indexFilename = BamIndexFactory::FindIndexFilename(reader->Filename(), preferredType);
204 // if no index file found (of any type)
205 if ( indexFilename.empty() ) {
206 const string message = string("could not find index file for:") + reader->Filename();
207 SetErrorString("BamRandomAccessController::LocateIndex", message);
211 // otherwise open & use index file that was found
212 return OpenIndex(indexFilename, reader);
215 bool BamRandomAccessController::OpenIndex(const string& indexFilename, BamReaderPrivate* reader) {
217 // attempt create new index of type based on filename
218 BamIndex* index = BamIndexFactory::CreateIndexFromFilename(indexFilename, reader);
220 const string message = string("could not open index file: ") + indexFilename;
221 SetErrorString("BamRandomAccessController::OpenIndex", message);
225 // attempt to load data from index file
226 if ( !index->Load(indexFilename) ) {
227 const string indexError = index->GetErrorString();
228 const string message = string("could not load index data from file: ") + indexFilename +
230 SetErrorString("BamRandomAccessController::OpenIndex", message);
234 // save new index & return success
239 bool BamRandomAccessController::RegionHasAlignments(void) const {
240 return m_hasAlignmentsInRegion;
243 void BamRandomAccessController::SetErrorString(const string& where, const string& what) {
244 m_errorString = where + ": " + what;
247 void BamRandomAccessController::SetIndex(BamIndex* index) {
253 bool BamRandomAccessController::SetRegion(const BamRegion& region, const int& referenceCount) {
258 // cannot jump when no index is available
260 SetErrorString("BamRandomAccessController", "cannot jump if no index data available");
264 // adjust region as necessary to reflect where data actually begins
265 AdjustRegion(referenceCount);
267 // if no data present, return true
268 // * Not an error, but future attempts to access alignments in this region will not return data
269 // Returning true is useful in a BamMultiReader setting where some BAM files may
270 // lack alignments in regions where other files still have data available.
271 if ( !m_hasAlignmentsInRegion )
274 // return success/failure of jump to specified region,
276 // * Index::Jump() is allowed to modify the m_hasAlignmentsInRegion flag
277 // This covers 'corner case' where a region is requested that lies beyond the last
278 // alignment on a reference. If this occurs, any subsequent calls to GetNextAlignment[Core]
279 // will not return data. BamMultiReader will still be able to successfully pull alignments
280 // from a region from other files even if this one has no data.
281 if ( !m_index->Jump(m_region, &m_hasAlignmentsInRegion) ) {
282 const string indexError = m_index->GetErrorString();
283 const string message = string("could not set region\n\t") + indexError;
284 SetErrorString("BamRandomAccessController::OpenIndex", message);