only writing histograms to files and this is separate from verbose

Author: andrewdavidsmith

src/bound_pop.cpp

@@ -47,8 +47,6 @@ using std::setprecision;
 using std::string;
 using std::vector;
 
-namespace fs = std::filesystem;
-
 static void
 report_bootstrapped_moments(const vector<double> &bootstrap_moments,
                             const MomentSequence &bootstrap_mom_seq,
@@ -112,7 +110,7 @@ bound_pop_main(const int argc, const char *argv[]) {
 Estimate a bound on the size of the underlying population based on
 counts of observed species in an initial sample.
 )";
-    string program_name = fs::path(argv[0]).filename();
+    string program_name = std::filesystem::path(argv[0]).filename();
     program_name += " " + string(argv[1]);
 
     /********** GET COMMAND LINE ARGUMENTS FOR BOUND_POP ***********/
@@ -239,13 +237,14 @@ counts of observed species in an initial sample.
            << "DISTINCT OBSERVATIONS  = " << distinct_obs << endl
            << "MAX COUNT              = " << counts_hist.size() - 1 << endl;
 
-      report_histogram(histogram_outfile, counts_hist);
-
       cerr << "OBSERVED MOMENTS" << endl;
       for (size_t i = 0; i < measure_moments.size(); i++)
         cerr << std::setprecision(16) << measure_moments[i] << endl;
     }
 
+    if (!histogram_outfile.empty())
+      report_histogram(histogram_outfile, counts_hist);
+
     if (QUICK_MODE) {
       if (measure_moments.size() < 2 * max_num_points)
         max_num_points = static_cast<size_t>(floor(measure_moments.size() / 2));

src/c_curve.cpp

@@ -46,8 +46,6 @@
 #include <string>
 #include <vector>
 
-namespace fs = std::filesystem;
-
 using std::accumulate;
 using std::array;
 using std::cbegin;
@@ -102,7 +100,7 @@ c_curve_main(const int argc, const char *argv[]) {
 Generate the complexity curve for data. This does not extrapolate, but
 instead resamples from the given data.
 )";
-    string program_name = fs::path(argv[0]).filename();
+    string program_name = std::filesystem::path(argv[0]).filename();
     program_name += " " + string(argv[1]);
 
     /********** GET COMMAND LINE ARGUMENTS  FOR C_CURVE ***********/
@@ -210,15 +208,17 @@ instead resamples from the given data.
       std::count_if(cbegin(counts_hist), cend(counts_hist),
                     [](const double x) { return x > 0.0; });
 
-    if (verbose) {
+    if (verbose)
       cerr << "TOTAL READS     = " << n_reads << endl
            << "COUNTS_SUM      = " << total_reads << endl
            << "DISTINCT READS  = " << distinct_reads << endl
            << "DISTINCT COUNTS = " << distinct_counts << endl
            << "MAX COUNT       = " << max_observed_count << endl
            << "COUNTS OF 1     = " << counts_hist[1] << endl;
+
+    if (!histogram_outfile.empty())
       report_histogram(histogram_outfile, counts_hist);
-    }
+
     const size_t upper_limit = n_reads;  // set upper limit equal to number of
                                          // molecules
 

src/common.hpp

@@ -24,7 +24,6 @@
 #include <cstddef>  // std::size_t
 #include <cstdint>  // std::uint64_t
 #include <fstream>
-#include <iostream>
 #include <ostream>
 #include <random>
 #include <stdexcept>
@@ -119,15 +118,12 @@ multinomial(std::mt19937 &gen, const std::vector<double> &mult_probs,
 template <typename H>
 void
 report_histogram(const std::string &outfile, const H &h) {
-  std::ofstream of;
-  if (!outfile.empty())
-    of.open(outfile);
-  std::ostream o(outfile.empty() ? std::cerr.rdbuf() : of.rdbuf());
-  o << "OBSERVED COUNTS (" << std::size(h) << ")" << std::endl;
+  std::ofstream out(outfile);
+  if (!out)
+    throw std::runtime_error("failed to open output file: " + outfile);
   for (auto i = 0u; i < std::size(h); ++i)
     if (h[i] > 0)
-      o << i << '\t' << static_cast<std::uint32_t>(h[i]) << std::endl;
-  o << std::endl;
+      out << i << '\t' << static_cast<std::uint32_t>(h[i]) << '\n';
 }
 
 #endif  // SRC_COMMON_HPP_

src/gc_extrap.cpp

@@ -42,8 +42,6 @@ using std::runtime_error;
 using std::string;
 using std::vector;
 
-namespace fs = std::filesystem;
-
 // ADS: functions same, header different (above and this one)
 static void
 write_predicted_coverage_curve(const string &outfile, const double c_level,
@@ -112,7 +110,7 @@ power-series expansion for the number of "unobserved" bases in the
 initial sample. The gc_extrap method is adapted to deal with
 individual nucleotides rather than distinct reads.
 )";
-    string program_name = fs::path(argv[0]).filename();
+    string program_name = std::filesystem::path(argv[0]).filename();
     program_name += " " + string(argv[1]);
 
     // ********* GET COMMAND LINE ARGUMENTS  FOR GC EXTRAP **********
@@ -219,7 +217,7 @@ individual nucleotides rather than distinct reads.
 
     orig_max_terms = min(orig_max_terms, first_zero - 1);
 
-    if (verbose) {
+    if (verbose)
       cerr << "TOTAL READS         = " << n_reads << endl
            << "BASE STEP SIZE      = " << base_step_size << endl
            << "BIN STEP SIZE       = " << bin_step_size << endl
@@ -230,8 +228,9 @@ individual nucleotides rather than distinct reads.
            << "TOTAL COVERED BASES = " << distinct_bins * bin_size << endl
            << "MAX COVERAGE COUNT  = " << max_observed_count << endl
            << "COUNTS OF 1         = " << coverage_hist[1] << endl;
+
+    if (!histogram_outfile.empty())
       report_histogram(histogram_outfile, coverage_hist);
-    }
 
     // catch if all reads are distinct
     if (orig_max_terms < MIN_REQUIRED_COUNTS)

src/lc_extrap.cpp

@@ -51,8 +51,6 @@ using std::uint32_t;
 using std::uint64_t;
 using std::vector;
 
-namespace fs = std::filesystem;
-
 int
 lc_extrap_main(const int argc, const char **argv) {
   try {
@@ -95,7 +93,7 @@ original goal of estimating the number of distinct reads that a
 sequencing library would yield upon deeper sequencing. This
 method has been used for many different purposes since then.
 )";
-    string program_name = fs::path(argv[0]).filename();
+    string program_name = std::filesystem::path(argv[0]).filename();
     program_name += " " + string(argv[1]);
 
     /********** GET COMMAND LINE ARGUMENTS  FOR LC EXTRAP ***********/
@@ -217,15 +215,16 @@ method has been used for many different purposes since then.
       std::count_if(cbegin(counts_hist), cend(counts_hist),
                     [](const double x) { return x > 0.0; });
 
-    if (verbose) {
+    if (verbose)
       cerr << "TOTAL READS     = " << n_reads << endl
            << "DISTINCT READS  = " << distinct_reads << endl
            << "DISTINCT COUNTS = " << distinct_counts << endl
            << "MAX COUNT       = " << max_observed_count << endl
            << "COUNTS OF 1     = " << counts_hist[1] << endl
            << "MAX TERMS       = " << orig_max_terms << endl;
+
+    if (!histogram_outfile.empty())
       report_histogram(histogram_outfile, counts_hist);
-    }
 
     // check to make sure library is not overly saturated
     const double two_fold_extrap = GoodToulmin2xExtrap(counts_hist);

src/pop_size.cpp

@@ -45,8 +45,6 @@ using std::string;
 using std::to_string;
 using std::vector;
 
-namespace fs = std::filesystem;
-
 int
 pop_size_main(const int argc, const char *argv[]) {
   try {
@@ -88,7 +86,7 @@ parameters assume that the initial sample represents at least 1e-9 of
 the population, which is sufficient for every example application we
 have seen.
 )";
-    string program_name = fs::path(argv[0]).filename();
+    string program_name = std::filesystem::path(argv[0]).filename();
     program_name += " " + string(argv[1]);
 
     /********** GET COMMAND LINE ARGUMENTS  FOR LC EXTRAP ***********/
@@ -213,15 +211,16 @@ have seen.
       std::count_if(begin(counts_hist), end(counts_hist),
                     [](const double x) { return x > 0.0; });
 
-    if (verbose) {
+    if (verbose)
       cerr << "TOTAL READS     = " << n_reads << endl
            << "DISTINCT READS  = " << distinct_reads << endl
            << "DISTINCT COUNTS = " << distinct_counts << endl
            << "MAX COUNT       = " << max_observed_count << endl
            << "COUNTS OF 1     = " << counts_hist[1] << endl
            << "MAX TERMS       = " << orig_max_terms << endl;
+
+    if (!histogram_outfile.empty())
       report_histogram(histogram_outfile, counts_hist);
-    }
 
     // check to make sure library is not overly saturated
     const double two_fold_extrap = GoodToulmin2xExtrap(counts_hist);