Changeset 3358

Timestamp:

Nov 23, 2014, 8:02:26 AM (8 years ago)

Author:

Peter

Message:

fix how to read from a hts_itr_t. refs #794

Location:

trunk/yat/omic

Files:

• BamFile.cc (modified) (1 diff)
• BamHeader.cc (modified) (2 diffs)
• BamReadIterator.cc (modified) (1 diff)

trunk/yat/omic/BamFile.cc

@@ -125 +125 @@
   bool InBamFile::read(BamRead& read, hts_itr_t* iter)
   {
-    assert(0 && "FIXME");
-    int result_ = -1;
+    int result = sam_itr_next(sf_, iter, read.bam_);
 #else
   bool InBamFile::read(BamRead& read, bam_iter_t iter)
   {
     assert(sf_->type & 1); // no random access from sam files
-    int result_ = bam_iter_read(sf_->x.bam, iter, read.bam_);
-#endif
-    if (result_<-1) {
+    int result = bam_iter_read(sf_->x.bam, iter, read.bam_);
+#endif
+    if (result<-1) {
       std::ostringstream ss;
       ss << "read(2): invalid bam file '" << filename() << "'";
-      if (result_ == -2)
+      if (result == -2)
         ss << " truncated file";
       else
-        ss << " error: " << result_;
+        ss << " error: " << result;
       throw utility::runtime_error(ss.str());
     }
-    return result_>=0;
+    return result>=0;
   }
 

trunk/yat/omic/BamHeader.cc

@@ -40 +40 @@
 
 
-  void BamHeader::parse_region(const std::string& reg, int& tid, int& begin,
+  void BamHeader::parse_region(const std::string& reg, int& id, int& begin,
                                int& end) const
   {
     assert(header_);
 #if YAT_HAVE_HTSLIB
-    assert(0 && "FIXME: parse_region not implement against htslib");
+    const char* b = reg.c_str();
+    const char* e = hts_parse_reg(b, &begin, &end);
+    // If begin > end suggests something is wrong and that is also how
+    // bam_parse_region (below in samtools impl) used to detect error
+    // and return non-zero.
+    if (begin<=end) {
+      std::string chr = reg.substr(0, e-b);
+      try {
+        id = tid(chr);
+        return;
+      }
+      catch (utility::runtime_error& e) {
+        ; // throw below instead
+      }
+    }
 #else
-    if (!bam_parse_region(header_, reg.c_str(), &tid, &begin, &end))
+    if (!bam_parse_region(header_, reg.c_str(), &id, &begin, &end))
       return;
 #endif
@@ -73 +87 @@
   {
 #if YAT_HAVE_HTSLIB
-    return bam_name2id(header_, name.c_str());
+    int res = bam_name2id(header_, name.c_str());
+    if (res==-1) {
+      std::string msg("invalid segment name:");
+      msg += name;
+      throw utility::runtime_error(msg);
+    }
+    return res;
 #else
     if (!header_->hash) {

trunk/yat/omic/BamReadIterator.cc

@@ -113 +113 @@
   {
 #if YAT_HAVE_HTSLIB
-    hts_itr_t* it = hts_itr_query(bf.index(), tid, begin, end, NULL);
+    hts_itr_t* it = sam_itr_queryi(bf.index(), tid, begin, end);
 #else
     bam_iter_t it = bam_iter_query(bf.index(), tid, begin, end);