multiple_variant_reader.h

Go to the documentation of this file.

#ifndef gamgee__multiple_variant_reader__guard
#define gamgee__multiple_variant_reader__guard

#include "htslib/vcf.h"

#include "variant_header.h"
#include "variant_header_merger.h"

#include "../exceptions.h"
#include "../utils/hts_memory.h"
#include "../utils/variant_utils.h"

namespace gamgee {

template<class ITERATOR>
class MultipleVariantReader {
 public:

  explicit MultipleVariantReader(const std::vector<std::string>& filenames, const bool validate_headers = true) :
    m_variant_files { },
    m_variant_headers { }
  {
    init_reader(filenames, validate_headers);
  }

  MultipleVariantReader(const std::vector<std::string>& filenames, const bool validate_headers,
                        const std::vector<std::string>& samples, const bool include = true) :
    m_variant_files { },
    m_variant_headers { }
  {
    init_reader(filenames, validate_headers);
    subset_variant_samples(m_variant_header_merger.get_raw_merged_header().get(), samples, include);
  }

  void init_reader(const std::vector<std::string>& filenames, const bool validate_headers) {
    m_variant_files.reserve(filenames.size());
    m_variant_headers.reserve(filenames.size());

    for (const auto& filename : filenames) {
      // TODO? check for maximum one stream
      auto* file_ptr = bcf_open(filename.empty() ? "-" : filename.c_str(), "r");
      if ( file_ptr == nullptr ) {
        throw FileOpenException{filename};
      }
      m_variant_files.push_back(std::move(utils::make_shared_hts_file(file_ptr)));

      auto* header_raw_ptr = bcf_hdr_read(file_ptr);
      if ( header_raw_ptr == nullptr ) {
        throw HeaderReadException{filename};
      }
      const auto& header_ptr = utils::make_shared_variant_header(header_raw_ptr);
      m_variant_headers.push_back(header_ptr);
 
      if (validate_headers && m_variant_header_merger.get_raw_merged_header().get() != nullptr)
  validate_header(header_ptr);  
      m_variant_header_merger.add_header(header_ptr); 
    }
  }

  MultipleVariantReader(MultipleVariantReader&& other) = default;
  MultipleVariantReader& operator=(MultipleVariantReader&& other) = default;

  MultipleVariantReader(const MultipleVariantReader&) = delete;
  MultipleVariantReader& operator=(const MultipleVariantReader& other) = delete;

  ITERATOR begin() const {
    return ITERATOR{m_variant_files, m_variant_headers};
  }

  ITERATOR end() const {
    return ITERATOR{};
  }

  const inline VariantHeader combined_header() const { return VariantHeader {m_variant_header_merger.get_raw_merged_header()}; }

  inline InputOrderedVariantHeaderMerger& get_variant_header_merger() { return m_variant_header_merger; }

  const std::vector<std::shared_ptr<bcf_hdr_t>>& get_input_vcf_headers() const { return m_variant_headers; }
  
 private:
  // TODO? only handles chromosome names, not lengths
  void validate_header(const std::shared_ptr<bcf_hdr_t>& other_header_ptr) {
    const auto& other_header = VariantHeader{other_header_ptr};
    if (combined_header().chromosomes() != other_header.chromosomes())
      throw HeaderCompatibilityException{"chromosomes in header files are inconsistent"};
  }

  std::vector<std::shared_ptr<htsFile>> m_variant_files;        
  std::vector<std::shared_ptr<bcf_hdr_t>> m_variant_headers;    
  InputOrderedVariantHeaderMerger m_variant_header_merger;      

};

}  // end namespace gamgee

#endif /* defined(gamgee__multiple_variant_reader__guard) */