source: trunk/yat/omic/BamFile.cc

Last change on this file was 3999, checked in by Peter, 13 months ago

update copyright years

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 5.1 KB
Line 
1// $Id: BamFile.cc 3999 2020-10-08 23:22:32Z peter $
2
3/*
4  Copyright (C) 2012, 2013, 2014, 2015, 2016, 2017, 2018, 2020 Peter Johansson
5
6  This file is part of the yat library, http://dev.thep.lu.se/yat
7
8  The yat library is free software; you can redistribute it and/or
9  modify it under the terms of the GNU General Public License as
10  published by the Free Software Foundation; either version 3 of the
11  License, or (at your option) any later version.
12
13  The yat library is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16  General Public License for more details.
17
18  You should have received a copy of the GNU General Public License
19  along with yat. If not, see <http://www.gnu.org/licenses/>.
20*/
21
22#include <config.h>
23
24#include "BamFile.h"
25#include "BamHeader.h"
26#include "BamReadIterator.h"
27
28#include "yat/utility/Exception.h"
29
30#include <htslib/sam.h>
31
32#include <cassert>
33#include <stdexcept>
34#include <string>
35#include <sstream>
36
37namespace theplu {
38namespace yat {
39namespace omic {
40
41  InBamFile::InBamFile(void)
42    : super_t(), index_(NULL)
43  {}
44
45
46  InBamFile::InBamFile(const std::string& fn)
47    : super_t(), index_(NULL)
48  {
49    open(fn);
50  }
51
52
53  InBamFile::~InBamFile(void)
54  {
55    hts_idx_destroy(index_);
56  }
57
58
59  void InBamFile::build_index(void) const
60  {
61    build_index_base();
62  }
63
64
65  uint64_t InBamFile::get_idx_stat(int tid, bool return_mapped) const
66  {
67    uint64_t mapped;
68    uint64_t unmapped;
69    if (tid >= header().n_targets())
70      throw utility::runtime_error("InBamFile::get_idx_stat");
71    hts_idx_get_stat(index(), tid, &mapped, &unmapped);
72    if (return_mapped)
73      return mapped;
74    return unmapped;
75  }
76
77
78  const BamHeader& InBamFile::header(void) const
79  {
80    return header_;
81  }
82
83
84  const InBamFile::index_type* InBamFile::index(void) const
85  {
86    if (!index_) {
87      index_ = hts_idx_load(filename().c_str(), HTS_FMT_BAI);
88      if (!index_) {
89        std::ostringstream ss;
90        ss << "cannot load bam index '" << filename() << ".bai'";
91        throw utility::runtime_error(ss.str());
92      }
93    }
94    return index_;
95  }
96
97
98  uint64_t InBamFile::n_mapped(int tid) const
99  {
100    return get_idx_stat(tid, true);
101  }
102
103
104  uint64_t InBamFile::n_unmapped(int tid) const
105  {
106    return get_idx_stat(tid, false);
107  }
108
109
110  uint64_t InBamFile::n_no_coordinate(void) const
111  {
112    return hts_idx_get_n_no_coor(index());
113  }
114
115
116  void InBamFile::open(const std::string& fn)
117  {
118    open_base(fn, "rb");
119    header_.header_ = sam_hdr_read(sf_);
120    if (header_.header_ == NULL) {
121      std::ostringstream ss;
122      ss << "failed to read header from '" << filename() << "'";
123      throw utility::runtime_error(ss.str());
124    }
125  }
126
127
128  bool InBamFile::read(BamRead& read)
129  {
130    int result_ = sam_read1(sf_, header_.header_, read.bam_);
131    if (result_<-1) {
132      std::ostringstream ss;
133      ss << "InBamFile::read(1): invalid bam file '" << filename() << "'";
134      if (result_ == -2)
135        ss << " truncated file";
136      else
137        ss << " error: " << result_;
138      throw utility::runtime_error(ss.str());
139    }
140    return result_>=0;
141  }
142
143
144  bool InBamFile::read(BamRead& read, hts_itr_t* iter)
145  {
146    int result = sam_itr_next(sf_, iter, read.bam_);
147    if (result<-1) {
148      std::ostringstream ss;
149      ss << "BamRead::read(2): invalid bam file '" << filename() << "'";
150      if (result == -2)
151        ss << " truncated file";
152      else
153        ss << " error: " << result;
154      throw utility::runtime_error(ss.str());
155    }
156    return result>=0;
157  }
158
159
160  OutBamFile::OutBamFile(void)
161    : super_t() {}
162
163
164  OutBamFile::OutBamFile(const std::string& fn, const BamHeader& header)
165  {
166    open(fn, header);
167  }
168
169
170  OutBamFile::OutBamFile(const std::string& fn, const BamHeader& header,
171                         unsigned int c)
172  {
173    open(fn, header, c);
174  }
175
176
177  void OutBamFile::build_index(void) const
178  {
179    if (is_open()) {
180      std::ostringstream msg;
181      msg << "failed building index for '" << filename() << "': file is open";
182      throw utility::runtime_error(msg.str());
183    }
184    build_index_base();
185  }
186
187
188  void OutBamFile::open(const std::string& fn, const BamHeader& h)
189  {
190    if (!h.header_)
191      throw utility::runtime_error("OutBamFile::open with null header");
192    open(fn, std::string("wb"), h);
193  }
194
195
196  void OutBamFile::open(const std::string& fn, const BamHeader& h,
197                        unsigned int compression)
198  {
199    std::string mode("wb");
200    if (compression > 9) {
201      std::stringstream oss;
202      oss << "OutBamFile::open( " << fn << ", <header>, " << compression
203          << "): incorrect compression level\n";
204      throw std::invalid_argument(oss.str());
205    }
206    if (compression)
207      mode.push_back('0'+compression);
208    else
209      mode = "wu";
210    open(fn, mode, h);
211  }
212
213
214  void OutBamFile::open(const std::string& fn, const std::string& mode,
215                        const BamHeader& h)
216  {
217    open_base(fn, mode);
218    if (sam_hdr_write(sf_, h.header_)) {
219      std::ostringstream ss;
220      ss << "failed open '" << fn << "'";
221      throw utility::runtime_error(ss.str());
222    }
223  }
224
225
226  void OutBamFile::write(const BamRead& read)
227  {
228    assert(read.bam_);
229    if (is_open() && bam_write1(sf_->fp.bgzf, read.bam_)>0)
230      return;
231    throw OutBamFile::error(read);
232  }
233
234
235  // OutBamFile::error
236
237  OutBamFile::error::error(const BamRead& r)
238    : IO_error("OutBamFile::write failed"), read_(r) {}
239
240
241  OutBamFile::error::~error(void) throw() {}
242
243
244  const BamRead& OutBamFile::error::read(void) const
245  { return read_; }
246
247}}}
Note: See TracBrowser for help on using the repository browser.