Changeset 3353

Timestamp:

Nov 22, 2014, 11:26:14 AM (8 years ago)

Author:

Peter

Message:

refs #794. First version that compiles against htslib. Alibi implementations in some cases though.

Location:

trunk/yat/omic

Files:

• BamFile.cc (modified) (8 diffs)
• BamFile.h (modified) (7 diffs)
• BamReadFilter.h (modified) (1 diff)
• BamReadIterator.cc (modified) (2 diffs)
• BamReadIterator.h (modified) (2 diffs)

trunk/yat/omic/BamFile.cc

@@ -53 +53 @@
   InBamFile::~InBamFile(void)
   {
+#if YAT_HAVE_HTSLIB
+    hts_idx_destroy(index_);
+#else
     bam_index_destroy(index_);
+#endif
   }
 
@@ -63 +67 @@
 
 
-  const bam_index_t* InBamFile::index(void) const
+  const InBamFile::index_type* InBamFile::index(void) const
   {
     if (!index_) {
+#if YAT_HAVE_HTSLIB
+      index_ = hts_idx_load(filename().c_str(), HTS_FMT_BAI);
+#else
       index_ = bam_index_load(filename().c_str());
+#endif
       if (!index_) {
         std::ostringstream ss;
@@ -79 +87 @@
   void InBamFile::open(const std::string& fn)
   {
+#if YAT_HAVE_HTSLIB
+    assert(0 && "FIXME");
+#else
     open_base(fn, "rb", NULL);
     header_.header_ = sf_->header;
+#endif
     if (header_.header_ == NULL) {
       std::ostringstream ss;
@@ -91 +103 @@
   bool InBamFile::read(BamRead& read)
   {
+#if YAT_HAVE_HTSLIB
+    int result_ = sam_read1(sf_, header_.header_, read.bam_);
+#else
     int result_ = samread(sf_, read.bam_);
+#endif
     if (result_<-1) {
       std::ostringstream ss;
@@ -105 +121 @@
 
 
+#if YAT_HAVE_HTSLIB
+  bool InBamFile::read(BamRead& read, hts_itr_t* iter)
+  {
+    assert(0 && "FIXME");
+    int result_ = -1;
+#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) {
       std::ostringstream ss;
@@ -141 +164 @@
   void OutBamFile::open(const std::string& fn, const BamHeader& h)
   {
+#if YAT_HAVE_HTSLIB
+    assert(0 && "FIXME");
+#else
     open_base(fn, "wb", h.header_);
+#endif
   }
 
@@ -155 +182 @@
       throw std::invalid_argument(oss.str());
     }
+#if YAT_HAVE_HTSLIB
+    assert(0 && "FIXME");
+#else
     mode.push_back('0'+compression);
     open_base(fn, mode, h.header_);
+#endif
   }
 
@@ -162 +193 @@
   void OutBamFile::write(const BamRead& read)
   {
+#if YAT_HAVE_HTSLIB
+    if (bam_write1(sf_->fp.bgzf, read.bam_)) {
+#else
     if (samwrite(sf_, read.bam_)<=0) {
+#endif
       throw OutBamFile::error(read);
     }

trunk/yat/omic/BamFile.h

@@ -80 +80 @@
        \see samopen
      */
+#ifndef YAT_HAVE_HTSLIB
     void open_base(const std::string& fn, const std::string& mode,
                    const void* aux);
+#endif
+
+#ifndef YAT_HAVE_HTSLIB
+    /// Old samtools 0.1.xx had a struct samfile_t, which in htslib has
+    /// changed to samFile
+    typedef samfile_t samFile;
+#endif
+
     /**
        bam file handler
      */
-    samfile_t* sf_;
+    samFile* sf_;
 
     /**
@@ -108 +117 @@
   public:
     /**
+       Struct describing bam index.
+
+       If compiled against htslib, \c hts_idx_t. If compiled against
+       libbam, \c bam_index_t.
+
+       \see index(void)
+     */
+#if YAT_HAVE_HTSLIB
+    typedef hts_idx_t index_type;
+#else
+    typedef bam_index_t index_type;
+#endif
+
+    /**
        \brief Default constructor
      */
@@ -139 +162 @@
        reads from stdin, this function throws.
      */
-    const bam_index_t* index(void) const;
+    const index_type* index(void) const;
 
     /**
@@ -162 +185 @@
        \return true on success
      */
+#if YAT_HAVE_HTSLIB
+    bool read(BamRead& read, hts_itr_t* iter);
+#else
     bool read(BamRead& read, bam_iter_t iter);
+#endif
   private:
     BamHeader header_;
-    mutable bam_index_t* index_;
+    mutable index_type* index_;
   };
 
@@ -280 +307 @@
   void BamFile<Derived>::close(void)
   {
+#if YAT_HAVE_HTSLIB
+    // FIXME, check that sam_close returns 0 (for success)
+    sam_close(sf_);
+#else
     samclose(sf_);
+#endif
     sf_ = NULL;
   }
@@ -292 +324 @@
 
 
+#ifndef YAT_HAVE_HTSLIB
   template<class Derived>
   void BamFile<Derived>::open_base(const std::string& fn,
@@ -306 +339 @@
     }
   }
+#endif
 
 }}}

trunk/yat/omic/BamReadFilter.h

@@ -27 +27 @@
 
 #include <functional>
+
+#ifndef BAM_DEF_MASK
+#define BAM_DEF_MASK (BAM_FUNMAP | BAM_FSECONDARY | BAM_FQCFAIL | BAM_FDUP)
+#endif
 
 namespace theplu {

trunk/yat/omic/BamReadIterator.cc

@@ -112 +112 @@
     : Actor(&bf)
   {
+#if YAT_HAVE_HTSLIB
+    hts_itr_t* it = hts_itr_query(bf.index(), tid, begin, end, NULL);
+#else
     bam_iter_t it = bam_iter_query(bf.index(), tid, begin, end);
+#endif
     if (!it) {
       std::stringstream ss;
@@ -120 +124 @@
     }
     // iter_ takes ownership
+#if YAT_HAVE_HTSLIB
+    iter_ = boost::shared_ptr<hts_itr_t>(it, IndexDestroyer());
+#else
     iter_ = boost::shared_ptr<__bam_iter_t>(it, IndexDestroyer());
+#endif
     assert(iter_.get());
   }

trunk/yat/omic/BamReadIterator.h

@@ -132 +132 @@
       void increment(void);
     private:
+#if YAT_HAVE_HTSLIB
+      boost::shared_ptr<hts_itr_t> iter_;
+#else
       boost::shared_ptr<__bam_iter_t> iter_;
+#endif
     };
 
@@ -139 +143 @@
     struct IndexDestroyer
     {
+#if YAT_HAVE_HTSLIB
+      void operator()(hts_itr_t* i) const
+      {
+        bam_itr_destroy(i);
+      }
+#else
       void operator()(bam_iter_t i) const
       {
         bam_iter_destroy(i);
       }
+#endif
     };
   };