fixing conflict with HEAD

Author: guilhermesena1

Makefile.am

@@ -25,7 +25,7 @@ lib_LIBRARIES = libsmithlab_cpp.a
 libsmithlab_cpp_a_SOURCES = GenomicRegion.cpp MappedRead.cpp		\
 OptionParser.cpp QualityScore.cpp bisulfite_utils.cpp			\
 chromosome_utils.cpp sim_utils.cpp smithlab_os.cpp smithlab_utils.cpp	\
-zlib_wrapper.cpp dna_four_bit.cpp cigar_utils.cpp
+zlib_wrapper.cpp dna_four_bit.cpp cigar_utils.cpp sam_record.cpp
 
 if ENABLE_HTS
 libsmithlab_cpp_a_SOURCES += htslib_wrapper_deprecated.cpp \

bisulfite_utils.hpp

@@ -29,30 +29,29 @@ namespace bsflags {
   // ADS: this is our addition to the SAM flags, using the "free" bits
   /* 4096 0x1000 read is A-rich
    */
-  static const uint16_t read_is_t_rich = 0x1000;
+  static const uint16_t read_is_a_rich = 0x1000;
 }
 
 inline bool
 is_t_rich(const sam_rec &sr) {
-  return samflags::check(sr.flags, bsflags::read_is_t_rich);
+  return !samflags::check(sr.flags, bsflags::read_is_a_rich);
 }
 
 inline bool
 is_a_rich(const sam_rec &sr) {
-  return !samflags::check(sr.flags, bsflags::read_is_t_rich);
+  return samflags::check(sr.flags, bsflags::read_is_a_rich);
 }
 
 inline void
 set_t_rich(sam_rec &sr) {
-  samflags::set(sr.flags, bsflags::read_is_t_rich);
+  samflags::unset(sr.flags, bsflags::read_is_a_rich);
 }
 
 inline void
 set_a_rich(sam_rec &sr) {
-  samflags::unset(sr.flags, bsflags::read_is_t_rich);
+  samflags::set(sr.flags, bsflags::read_is_a_rich);
 }
 
-
 void
 bisulfite_treatment(std::mt19937 &generator, std::string &seq,
                     double bs_rate = 1.0, double meth_rate = 0.0);

sam_record.cpp

@@ -33,11 +33,13 @@ using std::istream;
 using std::ostream;
 using std::begin;
 using std::end;
+using std::ostringstream;
+using std::to_string;
 
 // ADS: this is for debugging purposes
 string
 format_sam_flags(const uint16_t the_flags) {
-  std::ostringstream oss;
+  ostringstream oss;
   using samflags::check;
   oss << "read_paired: " << check(the_flags, samflags::read_paired) << '\n'
       << "read_pair_mapped: "
@@ -72,23 +74,36 @@ format_sam_flags(const uint16_t the_flags) {
 //   return regex_match(qual, regex("[!-~]+"));
 // }
 
+size_t
+sam_rec::estimate_line_size() const {
+  static const size_t all_field_estimates = 100;
+  return qname.size() + rname.size() + qual.size() + all_field_estimates;
+}
+
+string
+sam_rec::tostring() const {
+  string out;
+  out.reserve(estimate_line_size());
+  out.append(qname + "\t" +
+               to_string(flags) + "\t" +
+               rname + "\t" +
+               to_string(pos) + "\t" +
+               to_string(static_cast<unsigned>(mapq)) + "\t" +
+               cigar + "\t" +
+               rnext + "\t" +
+               to_string(pnext) + "\t" +
+               to_string(tlen) + "\t" +
+               seq + "\t" +
+               qual);
+  for (auto it(begin(tags)); it != end(tags); ++it)
+    out.append("\t" + *it);
+
+  return out;
+}
 
 ostream &
 operator<<(std::ostream &the_stream, const sam_rec &r) {
-  the_stream << r.qname << '\t'
-             << r.flags << '\t'
-             << r.rname << '\t'
-             << r.pos << '\t'
-             << static_cast<unsigned>(r.mapq) << '\t'
-             << r.cigar << '\t'
-             << r.rnext << '\t'
-             << r.pnext << '\t'
-             << r.tlen << '\t'
-             << r.seq << '\t'
-             << r.qual;
-
-  for (auto it(begin(r.tags)); it != end(r.tags); ++it)
-      the_stream << '\t' << *it;
+  the_stream << r.tostring();
   return the_stream;
 }
 

sam_record.hpp

@@ -22,6 +22,8 @@
 #include <string>
 #include <vector>
 #include <iostream>
+#include <sstream>
+#include <iterator>
 
 // from 30 April 2020 SAM documentation
 // 1    0x1   template having multiple segments in sequencing
@@ -117,6 +119,8 @@ class sam_rec {
     seq(_seq),
     qual(_qual) {}
   void add_tag(const std::string &the_tag) {tags.push_back(the_tag);}
+  size_t estimate_line_size() const;
+  std::string tostring() const;
 };
 
 inline bool
@@ -144,4 +148,36 @@ void
 inflate_with_cigar(const sam_rec &sr, std::string &to_inflate,
                    const char inflation_symbol = 'N');
 
+template<typename T>
+static void
+write_sam_header(const std::vector<std::string> &chrom_names,
+                 const std::vector<T> &chrom_starts,
+                 const std::string program_name,
+                 const std::string program_version,
+                 const int argc, const char **argv,
+                 std::ostream &out) {
+  static const std::string SAM_VERSION = "1.0";
+
+  // sam version
+  out <<"@HD" << '\t' << "VN:" << SAM_VERSION << '\n'; // sam version
+
+  // chromosome sizes
+  const size_t n_chroms = chrom_names.size() - 1;
+  for (size_t i = 1; i < n_chroms; ++i) {
+    out << "@SQ" << '\t'
+        << "SN:" << chrom_names[i] << '\t'
+        << "LN:" << chrom_starts[i+1] - chrom_starts[i] << '\n';
+  }
+
+  // program details
+  out << "@PG" << '\t'
+      << "ID:" << program_name << '\t'
+      << "VN:" << program_version << '\t';
+
+  // how the program was run
+  std::ostringstream the_command;
+  copy(argv, argv + argc, std::ostream_iterator<const char*>(the_command, " "));
+  out << "CL:\"" << the_command.str() << "\"" << std::endl;
+}
+
 #endif