src/api/internal/BamRandomAccessController_p.cpp

   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 // **************************************************************************
   9
  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;
  17
  18 #include <cassert>
  19 #include <sstream>
  20 using namespace std;
  21
  22 BamRandomAccessController::BamRandomAccessController(void)
  23     : m_index(0)
  24     , m_indexCacheMode(BamIndex::LimitedIndexCaching)
  25     , m_hasAlignmentsInRegion(true)
  26 { }
  27
  28 BamRandomAccessController::~BamRandomAccessController(void) {
  29     Close();
  30 }
  31
  32 void BamRandomAccessController::AdjustRegion(const int& referenceCount) {
  33
  34     // skip if no index available
  35     if ( m_index == 0 )
  36         return;
  37
  38     // see if any references in region have alignments
  39     m_hasAlignmentsInRegion = false;
  40     int currentId = m_region.LeftRefID;
  41     const int rightBoundRefId = ( m_region.isRightBoundSpecified() ? m_region.RightRefID : referenceCount - 1 );
  42     while ( currentId <= rightBoundRefId ) {
  43         m_hasAlignmentsInRegion = m_index->HasAlignments(currentId);
  44         if ( m_hasAlignmentsInRegion ) break;
  45         ++currentId;
  46     }
  47
  48     // if no data found on any reference in region
  49     if ( !m_hasAlignmentsInRegion )
  50         return;
  51
  52     // if left bound of desired region had no data, use first reference that had data
  53     // otherwise, leave requested region as-is
  54     if ( currentId != m_region.LeftRefID ) {
  55         m_region.LeftRefID = currentId;
  56         m_region.LeftPosition = 0;
  57     }
  58 }
  59
  60 // returns alignments' "RegionState": { Before|Overlaps|After } current region
  61 BamRandomAccessController::RegionState
  62 BamRandomAccessController::AlignmentState(const BamAlignment& alignment) const {
  63
  64     // if region has no left bound at all
  65     if ( !m_region.isLeftBoundSpecified() )
  66         return OverlapsRegion;
  67
  68     // handle unmapped reads - return AFTER region to halt processing
  69     if ( alignment.RefID == -1 )
  70         return AfterRegion;
  71
  72     // if alignment is on any reference before left bound reference
  73     if ( alignment.RefID < m_region.LeftRefID )
  74         return BeforeRegion;
  75
  76     // if alignment is on left bound reference
  77     else if ( alignment.RefID == m_region.LeftRefID ) {
  78
  79         // if alignment starts at or after left bound position
  80         if ( alignment.Position >= m_region.LeftPosition) {
  81
  82             if ( m_region.isRightBoundSpecified() &&             // right bound is specified AND
  83                  m_region.LeftRefID == m_region.RightRefID &&    // left & right bounds on same reference AND
  84                  alignment.Position >= m_region.RightPosition )  // alignment starts on or after right bound position
  85                 return AfterRegion;
  86
  87             // otherwise, alignment overlaps region
  88             else return OverlapsRegion;
  89         }
  90
  91         // alignment starts before left bound position
  92         else {
  93
  94             // if alignment overlaps left bound position
  95             if ( alignment.GetEndPosition() > m_region.LeftPosition )
  96                 return OverlapsRegion;
  97             else
  98                 return BeforeRegion;
  99         }
 100     }
 101
 102     // otherwise alignment is on a reference after left bound reference
 103     else {
 104
 105         // if region has a right bound
 106         if ( m_region.isRightBoundSpecified() ) {
 107
 108             // alignment is on any reference between boundaries
 109             if ( alignment.RefID < m_region.RightRefID )
 110                 return OverlapsRegion;
 111
 112             // alignment is on any reference after right boundary
 113             else if ( alignment.RefID > m_region.RightRefID )
 114                 return AfterRegion;
 115
 116             // alignment is on right bound reference
 117             else {
 118
 119                 // if alignment starts before right bound position
 120                 if ( alignment.Position < m_region.RightPosition )
 121                     return OverlapsRegion;
 122                 else
 123                     return AfterRegion;
 124             }
 125         }
 126
 127         // otherwise, alignment starts after left bound and there is no right bound given
 128         else return OverlapsRegion;
 129     }
 130 }
 131
 132 void BamRandomAccessController::Close(void) {
 133     ClearIndex();
 134     ClearRegion();
 135 }
 136
 137 void BamRandomAccessController::ClearIndex(void) {
 138     if ( m_index ) {
 139         delete m_index;
 140         m_index = 0;
 141     }
 142 }
 143
 144 void BamRandomAccessController::ClearRegion(void) {
 145     m_region.clear();
 146     m_hasAlignmentsInRegion = true;
 147 }
 148
 149 bool BamRandomAccessController::CreateIndex(BamReaderPrivate* reader,
 150                                             const BamIndex::IndexType& type)
 151 {
 152     // skip if reader is invalid
 153     assert(reader);
 154     if ( !reader->IsOpen() ) {
 155         SetErrorString("BamRandomAccessController::CreateIndex",
 156                        "cannot create index for unopened reader");
 157         return false;
 158     }
 159
 160     // create new index of requested type
 161     BamIndex* newIndex = BamIndexFactory::CreateIndexOfType(type, reader);
 162     if ( newIndex == 0 ) {
 163         stringstream s("");
 164         s << "could not create index of type: " << type;
 165         SetErrorString("BamRandomAccessController::CreateIndex", s.str());
 166         return false;
 167     }
 168
 169     // attempt to build index from current BamReader file
 170     if ( !newIndex->Create() ) {
 171         const string indexError = newIndex->GetErrorString();
 172         const string message = "could not create index: \n\t" + indexError;
 173         SetErrorString("BamRandomAccessController::CreateIndex", message);
 174         return false;
 175     }
 176
 177     // save new index
 178     SetIndex(newIndex);
 179
 180     // set new index's cache mode & return success
 181     newIndex->SetCacheMode(m_indexCacheMode);
 182     return true;
 183 }
 184
 185 string BamRandomAccessController::GetErrorString(void) const {
 186     return m_errorString;
 187 }
 188
 189 bool BamRandomAccessController::HasIndex(void) const {
 190     return ( m_index != 0 );
 191 }
 192
 193 bool BamRandomAccessController::HasRegion(void) const  {
 194     return ( !m_region.isNull() );
 195 }
 196
 197 bool BamRandomAccessController::IndexHasAlignmentsForReference(const int& refId) {
 198     return m_index->HasAlignments(refId);
 199 }
 200
 201 bool BamRandomAccessController::LocateIndex(BamReaderPrivate* reader,
 202                                             const BamIndex::IndexType& preferredType)
 203 {
 204     // look up index filename, deferring to preferredType if possible
 205     assert(reader);
 206     const string& indexFilename = BamIndexFactory::FindIndexFilename(reader->Filename(), preferredType);
 207
 208     // if no index file found (of any type)
 209     if ( indexFilename.empty() ) {
 210         const string message = string("could not find index file for:") + reader->Filename();
 211         SetErrorString("BamRandomAccessController::LocateIndex", message);
 212         return false;
 213     }
 214
 215     // otherwise open & use index file that was found
 216     return OpenIndex(indexFilename, reader);
 217 }
 218
 219 bool BamRandomAccessController::OpenIndex(const string& indexFilename, BamReaderPrivate* reader) {
 220
 221     // attempt create new index of type based on filename
 222     BamIndex* index = BamIndexFactory::CreateIndexFromFilename(indexFilename, reader);
 223     if ( index == 0 ) {
 224         const string message = string("could not open index file: ") + indexFilename;
 225         SetErrorString("BamRandomAccessController::OpenIndex", message);
 226         return false;
 227     }
 228
 229     // set cache mode
 230     index->SetCacheMode(m_indexCacheMode);
 231
 232     // attempt to load data from index file
 233     if ( !index->Load(indexFilename) ) {
 234         const string indexError = index->GetErrorString();
 235         const string message = string("could not load index data from file: ") + indexFilename +
 236                                "\n\t" + indexError;
 237         SetErrorString("BamRandomAccessController::OpenIndex", message);
 238         return false;
 239     }
 240
 241     // save new index & return success
 242     SetIndex(index);
 243     return true;
 244 }
 245
 246 bool BamRandomAccessController::RegionHasAlignments(void) const {
 247     return m_hasAlignmentsInRegion;
 248 }
 249
 250 void BamRandomAccessController::SetErrorString(const string& where, const string& what) {
 251     m_errorString = where + ": " + what;
 252 }
 253
 254 void BamRandomAccessController::SetIndex(BamIndex* index) {
 255     if ( m_index )
 256         ClearIndex();
 257     m_index = index;
 258 }
 259
 260 void BamRandomAccessController::SetIndexCacheMode(const BamIndex::IndexCacheMode& mode) {
 261     m_indexCacheMode = mode;
 262     if ( m_index )
 263         m_index->SetCacheMode(mode);
 264 }
 265
 266 bool BamRandomAccessController::SetRegion(const BamRegion& region, const int& referenceCount) {
 267
 268     // store region
 269     m_region = region;
 270
 271     // cannot jump when no index is available
 272     if ( !HasIndex() ) {
 273         SetErrorString("BamRandomAccessController", "cannot jump if no index data available");
 274         return false;
 275     }
 276
 277     // adjust region as necessary to reflect where data actually begins
 278     AdjustRegion(referenceCount);
 279
 280     // if no data present, return true
 281     //   * Not an error, but future attempts to access alignments in this region will not return data
 282     //     Returning true is useful in a BamMultiReader setting where some BAM files may
 283     //     lack alignments in regions where other files still have data available.
 284     if ( !m_hasAlignmentsInRegion )
 285         return true;
 286
 287     // return success/failure of jump to specified region,
 288     //
 289     //  * Index::Jump() is allowed to modify the m_hasAlignmentsInRegion flag
 290     //    This covers 'corner case' where a region is requested that lies beyond the last
 291     //    alignment on a reference. If this occurs, any subsequent calls to GetNextAlignment[Core]
 292     //    will not return data. BamMultiReader will still be able to successfully pull alignments
 293     //    from a region from other files even if this one has no data.
 294     if ( !m_index->Jump(m_region, &m_hasAlignmentsInRegion) ) {
 295         const string indexError = m_index->GetErrorString();
 296         const string message = string("could not set region\n\t") + indexError;
 297         SetErrorString("BamRandomAccessController::OpenIndex", message);
 298         return false;
 299     }
 300     else
 301         return true;
 302 }