diff --git a/CMakeLists.txt b/CMakeLists.txt
index 6604e83..827ae53 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -79,6 +79,7 @@ endif()
 
 add_library(adapter_set OBJECT src/adapter_set.cpp)
 add_library(falco_utils OBJECT src/falco_utils.cpp)
+add_library(falco_gc_summary OBJECT src/falco_gc_summary.cpp)
 add_library(base_groups OBJECT src/base_groups.cpp)
 add_library(bgzf_block OBJECT src/bgzf_block.cpp)
 add_library(bgzf_reader OBJECT src/bgzf_reader.cpp)
@@ -92,6 +93,7 @@ add_library(bam_file OBJECT src/bam_file.cpp)
 add_library(samrec OBJECT src/samrec.cpp)
 add_library(sam_file OBJECT src/sam_file.cpp)
 add_library(sam_stdin OBJECT src/sam_stdin.cpp)
+add_library(file_info OBJECT src/file_info.cpp)
 add_library(falco_file_format OBJECT src/falco_file_format.cpp)
 add_library(quality_score OBJECT src/quality_score.cpp)
 add_library(tile_processor OBJECT src/tile_processor.cpp)
@@ -118,6 +120,7 @@ target_link_libraries(fastq_bgzf_file PRIVATE
 target_link_libraries(bam_file PRIVATE HTSLIB::HTSLIB)
 target_link_libraries(falco_file_format PRIVATE HTSLIB::HTSLIB)
 target_link_libraries(tile_processor PRIVATE HTSLIB::HTSLIB)
+target_link_libraries(quality_score PRIVATE HTSLIB::HTSLIB)
 target_link_libraries(fastq_gz_file PRIVATE HTSLIB::HTSLIB)
 if(USE_ISAL)
   target_link_libraries(fastq_gz_file PUBLIC ISAL::ISAL)
@@ -142,8 +145,10 @@ target_link_libraries(falco PRIVATE
   sam_stdin
   samrec
   bam_header
+  file_info
   falco_file_format
   falco_utils
+  falco_gc_summary
   falco_grade
   falco_analyzer
   results_summary
diff --git a/src/bam_file.cpp b/src/bam_file.cpp
index 27d0acc..90533a8 100644
--- a/src/bam_file.cpp
+++ b/src/bam_file.cpp
@@ -123,15 +123,14 @@ bam_file::load_next(const std::int32_t file_id,
     input_last = n_to_keep;
   }
   br.release();  // use like monotonic_buffer_resource
-  // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic)
-  auto in_itr = std::data(input_buffer) + input_last;
+  auto in_itr = std::begin(input_buffer) + input_last;
   while (br.task_ready() && has_in()) {
-    auto task = br.get_decomp_task(in_itr);
-    in_itr += task.size;  // NOLINT(*-pro-bounds-pointer-arithmetic)
+    auto task = br.get_decomp_task(std::to_address(in_itr));
+    in_itr += task.size;
     ++n_tasks;
     tq.push(file_id, std::move(task));
   }
-  input_last = std::distance(std::data(input_buffer), in_itr);
+  input_last = std::distance(std::begin(input_buffer), in_itr);
 }
 
 [[nodiscard]] static inline auto
@@ -143,10 +142,11 @@ partition(auto itr,
           std::atomic_int32_t &n_tasks) {
   // ADS: this isn't working as desired: the end position of each part should be
   // the first record end past the 'end_itr' below unless end_itr == end
-  auto dist = std::distance(itr, end);
-  const auto chunk_size = (dist + n_chunks - 1) / n_chunks;
+  static constexpr std::int64_t min_chunk_size = 64 * 1024;
+  auto chunk_size = (std::distance(itr, end) + n_chunks - 1) / n_chunks;
+  chunk_size = std::max(min_chunk_size, chunk_size);
   while (itr != end) {
-    dist = std::distance(itr, end);
+    const auto dist = std::distance(itr, end);
     const auto end_itr = itr + (dist < chunk_size ? dist : chunk_size);
     // ADS: find_end_pos doesn't find end pos of a record, but of a range, so
     // includes multiple records
diff --git a/src/bam_file.hpp b/src/bam_file.hpp
index 00a7f27..97a2473 100644
--- a/src/bam_file.hpp
+++ b/src/bam_file.hpp
@@ -19,7 +19,7 @@
 struct task_queue;
 
 class bam_file {
-  static constexpr std::int64_t min_buf_size = 64 * 1024;
+  static constexpr std::int64_t min_buf_size = 256 * 1024;
   std::vector<char> input_buffer;
   std::vector<char> output_buffer;
   std::int64_t input_last{};
@@ -104,7 +104,7 @@ class bam_file {
                      std::atomic_int32_t &n_tasks) -> void;
 
   auto
-  make_tasks(const std::int64_t n_chunks,
+  make_tasks(const std::int64_t n_threads,
              const std::int32_t file_id,
              task_queue &tq,
              std::atomic_int32_t &n_tasks) -> void;
diff --git a/src/bamrec.cpp b/src/bamrec.cpp
index ed9eb6a..be29899 100644
--- a/src/bamrec.cpp
+++ b/src/bamrec.cpp
@@ -86,8 +86,7 @@ bamrec::get_next(bamrec::pos_t &itr, const bamrec::pos_t end, bamrec &rec)
   else
     assign_sequence(seq_in, core.l_seq, out_itr);
   out_itr += rec.seq_len;  // increment data cursor to qual
-  const auto has_qual = (itr[core.qual_offset()] != qual_missing_code);
-  if (has_qual) {
+  if (itr[core.qual_offset()] != qual_missing_code) {
     const auto qual_itr = itr + core.qual_offset();
     const auto qual_end = qual_itr + rec.seq_len;
     if (core.bam_is_rev())
diff --git a/src/bamrec.hpp b/src/bamrec.hpp
index a2fd618..0203a3e 100644
--- a/src/bamrec.hpp
+++ b/src/bamrec.hpp
@@ -80,6 +80,7 @@ class bamrec {
   friend constexpr auto get_seq(const bamrec &);
   friend constexpr auto get_seq_end(const bamrec &);
   friend constexpr auto get_seq_size(const bamrec &);
+  friend constexpr auto has_qual(const bamrec &);
   friend constexpr auto get_qual(const bamrec &);
   friend constexpr auto get_qual_end(const bamrec &);
   friend constexpr auto get_qual_size(const bamrec &);
@@ -102,6 +103,11 @@ class bamrec {
 
   [[nodiscard]] static auto
   find_end_pos(pos_t itr, const pos_t end) -> pos_t;
+
+  [[nodiscard]] static consteval auto
+  get_qual_missing_code() -> char {
+    return qual_missing_code;
+  }
 };
 
 [[nodiscard]] inline constexpr auto
@@ -135,13 +141,18 @@ get_qual(const bamrec &rec) {
 }
 
 [[nodiscard]] inline constexpr auto
-get_qual_end(const bamrec &rec) {
-  return get_qual(rec) + rec.seq_len;
+has_qual(const bamrec &rec) {
+  return *get_qual(rec) != bamrec::qual_missing_code;
 }
 
 [[nodiscard]] inline constexpr auto
 get_qual_size(const bamrec &rec) {
-  return get_seq_size(rec);
+  return has_qual(rec) ? rec.seq_len : 0U;
+}
+
+[[nodiscard]] inline constexpr auto
+get_qual_end(const bamrec &rec) {
+  return get_qual(rec) + get_qual_size(rec);
 }
 
 struct bam_task_t {
diff --git a/src/bgzf_block.hpp b/src/bgzf_block.hpp
index ef4ff61..1dc0157 100644
--- a/src/bgzf_block.hpp
+++ b/src/bgzf_block.hpp
@@ -8,15 +8,19 @@
 static constexpr auto max_bgzf_block_size = 65536;
 
 struct bgzf_block_t {
+  using iterator = char *;
+
   // ADS: because out_itr and in_itr reference separate blocks, either we need
   // to track one of the sizes or we would have to add an end iterator for one
   // of them.
   std::int32_t size{};
-  char *out_itr{};
-  char *in_itr{};
-
-  bgzf_block_t(const std::int32_t size, char *out_itr, char *in_itr) noexcept :
+  iterator out_itr{};
+  iterator in_itr{};
+  bgzf_block_t(const std::int32_t size,
+               const iterator out_itr,
+               const iterator in_itr) noexcept :
     size{size}, out_itr{out_itr}, in_itr{in_itr} {}
+
   // clang-format off
   bgzf_block_t(const bgzf_block_t &src) = delete;
   auto operator=(const bgzf_block_t &src) -> bgzf_block_t & = delete;
@@ -24,7 +28,6 @@ struct bgzf_block_t {
   ~bgzf_block_t() = default;
   bgzf_block_t(bgzf_block_t &&src) noexcept = default;
   auto operator=(bgzf_block_t &&src) noexcept -> bgzf_block_t & = default;
-  [[nodiscard]] auto data() const -> const char * { return in_itr; }
   operator bool() const { return size > 0; }
   // clang-format on
 
diff --git a/src/bgzf_reader.cpp b/src/bgzf_reader.cpp
index fcd3ceb..efe6fd6 100644
--- a/src/bgzf_reader.cpp
+++ b/src/bgzf_reader.cpp
@@ -24,29 +24,30 @@ bgzf_reader::bgzf_reader(const std::string &filename,
   filesize{std::filesystem::file_size(filename)},             //
   inbuf(std::make_unique_for_overwrite<char[]>(inbuf_size)),  // NOLINT
   outbuf(std::make_unique_for_overwrite<char[]>(buf_size)),   // NOLINT
-  next_in{inbuf.get()},                                       //
-  end_in{inbuf.get()},                                        //
-  next_out{outbuf.get()},                                     //
+  next_in_itr{inbuf.get()},                                   //
+  end_in_itr{inbuf.get()},                                    //
+  next_out_itr{outbuf.get()},                                 //
   // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic)
-  end_out{outbuf.get() + buf_size}  //
+  end_out_itr{outbuf.get() + buf_size}  //
 {}
 
 [[nodiscard]] auto
 bgzf_reader::read_data() -> bool {
   if (at_eof())
     return false;
-  const auto unused_in = std::distance(next_in, end_in);
-  std::copy_n(next_in, unused_in, inbuf.get());
-  next_in = inbuf.get();
+  const auto unused_in = std::distance(next_in_itr, end_in_itr);
+  std::copy_n(next_in_itr, unused_in, inbuf.get());
+  next_in_itr = inbuf.get();
   // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic)
-  end_in = next_in + unused_in;
+  end_in_itr = next_in_itr + unused_in;
   const auto avail_in = inbuf_size - unused_in;
-  const auto n_bytes = std::fread(end_in, 1, avail_in, fp.get());
+  const auto n_bytes = std::fread(end_in_itr, 1, avail_in, fp.get());
   if (std::ferror(fp.get()))
     throw std::system_error(std::make_error_code(std::errc(errno)),
                             "failed reading input");
   // will usually be end of inbuf
-  end_in += n_bytes;  // NOLINT(cppcoreguidelines-pro-bounds-pointer-arithmetic)
+  // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic)
+  end_in_itr += n_bytes;
   return n_bytes > 0;
 }
 
@@ -109,8 +110,8 @@ get_isize(const auto data, const auto data_size) {
 }
 
 static inline auto
-assign(gzip_header &hdr,  // cppcheck-suppress constParameterReference
-       const auto data) -> void {
+assign_gzip_hdr(gzip_header &hdr,  // cppcheck-suppress constParameterReference
+                const auto data) -> void {
   std::memcpy(std::data(hdr.data), std::to_address(data), gzip_header_size);
 }
 
@@ -123,22 +124,24 @@ get_gzip_body_size(const gzip_header &gh) -> std::uint32_t {
 
 [[nodiscard]] auto
 bgzf_reader::get_decomp_task(char *out_itr) -> bgzf_block_t {
-  if (std::distance(next_out, end_out) < max_bgzf_block_size)
+  if (std::distance(next_out_itr, end_out_itr) < max_bgzf_block_size)
     return bgzf_block_t{};
-  if (std::distance(next_in, end_in) < gzip_header_size && !read_data())
+  if (std::distance(next_in_itr, end_in_itr) < gzip_header_size && !read_data())
     return bgzf_block_t{};
   gzip_header gh;
-  assign(gh, next_in);
+  assign_gzip_hdr(gh, next_in_itr);
   assert(gh.check_magic());
   // NOLINTBEGIN(cppcoreguidelines-pro-bounds-pointer-arithmetic)
-  next_in += gzip_header_size;
+  next_in_itr += gzip_header_size;
   const auto body_size = get_gzip_body_size(gh);
-  if (std::distance(next_in, end_in) < body_size && !read_data())
+  if (std::distance(next_in_itr, end_in_itr) < body_size && !read_data())
     return bgzf_block_t{};
-  std::memcpy(next_out, next_in, body_size);
-  bgzf_block_t task(get_isize(next_in, body_size), out_itr, next_out);
-  next_in += body_size;
-  next_out += max_bgzf_block_size;
+  std::copy_n(next_in_itr, body_size, next_out_itr);
+  const auto isize = get_isize(next_in_itr, body_size);
+  assert(isize >= 0);
+  bgzf_block_t task(isize, out_itr, next_out_itr);
+  next_in_itr += body_size;
+  next_out_itr += max_bgzf_block_size;
   // NOLINTEND(cppcoreguidelines-pro-bounds-pointer-arithmetic)
   return task;
 }
diff --git a/src/bgzf_reader.hpp b/src/bgzf_reader.hpp
index dc148a2..276a35c 100644
--- a/src/bgzf_reader.hpp
+++ b/src/bgzf_reader.hpp
@@ -13,6 +13,9 @@
 
 // reads data and provides serialized compressed chunks for deflation
 class bgzf_reader {
+public:
+  using iterator = char *;
+
 private:
   static constexpr auto inbuf_size = 16 * max_bgzf_block_size;
 
@@ -20,10 +23,10 @@ class bgzf_reader {
   std::uint64_t filesize{};
   std::unique_ptr<char[]> inbuf;   // NOLINT(cppcoreguidelines-avoid-c-arrays)
   std::unique_ptr<char[]> outbuf;  // NOLINT(cppcoreguidelines-avoid-c-arrays)
-  char *next_in{};
-  char *end_in{};
-  char *next_out{};
-  char *end_out{};
+  iterator next_in_itr{};
+  iterator end_in_itr{};
+  iterator next_out_itr{};
+  iterator end_out_itr{};
 
 public:
   bgzf_reader(const std::string &filename, const std::int64_t buf_size);
@@ -31,7 +34,7 @@ class bgzf_reader {
   operator bool() const { return can_produce_data(); }
 
   [[nodiscard]] auto
-  get_decomp_task(char *out_itr) -> bgzf_block_t;
+  get_decomp_task(iterator out_itr) -> bgzf_block_t;
 
   [[nodiscard]] auto
   read_data() -> bool;
@@ -39,7 +42,7 @@ class bgzf_reader {
   auto
   release() {
     // ADS: intended analogous to monotonic_buffer_resource release
-    next_out = outbuf.get();
+    next_out_itr = outbuf.get();
   }
 
   [[nodiscard]] auto
@@ -64,12 +67,12 @@ class bgzf_reader {
 
   [[nodiscard]] auto
   can_produce_data() const -> bool {
-    return !at_eof() || next_in != end_in;
+    return !at_eof() || next_in_itr != end_in_itr;
   }
 
   [[nodiscard]] auto
   has_out() const -> bool {
-    return std::distance(next_out, end_out) >= max_bgzf_block_size;
+    return std::distance(next_out_itr, end_out_itr) >= max_bgzf_block_size;
   }
 };
 
diff --git a/src/falco.cpp b/src/falco.cpp
index e0a4c3f..d97cb05 100644
--- a/src/falco.cpp
+++ b/src/falco.cpp
@@ -66,6 +66,7 @@ Use these as templates. Copy and modify them to customize your analysis.
 #include "file_info.hpp"
 #include "get_binary_dir.hpp"
 #include "original_duplicates.hpp"
+#include "quality_score.hpp"
 #include "reads_file.hpp"  // IWYU pragma: keep
 #include "results_summary.hpp"
 #include "run_mode.hpp"
@@ -192,16 +193,20 @@ get_file_info_stdin(const std::vector<std::string> &names,
     ft_tile.first == falco::file_format::fastq
       ? estimate_n_reads_fastq_stdin(names.front())
       : estimate_n_reads_sam_stdin(names.front());
-  file_info info;
-  info.name = names.front();
-  info.format = ft_tile.first;
-  info.description = std::format("{} from standard input", info.format);
-  info.size = 0;
-  info.has_tiles = (ft_tile.second != 0);
-  info.tile_id_position = ft_tile.second;
-  info.n_reads_est = n_reads_est;
-  info.read_len_est = read_len_est;
-  return std::vector{1, info};
+  // clang-format off
+  return std::vector{1, file_info{
+    .name = names.front(),
+    .format = ft_tile.first,
+    .description = std::format("{} from standard input", ft_tile.first),
+    .size = 0,
+    .n_reads_est = n_reads_est,
+    .read_len_est = read_len_est,
+    .has_quals = true,
+    .encoding = falco::encoding::unknown,
+    .has_tiles = (ft_tile.second != 0),
+    .tile_id_position = ft_tile.second,
+    }};
+  // clang-format on
 }
 
 [[nodiscard]] static auto
@@ -211,6 +216,8 @@ get_file_info(const auto &infiles) {
     const auto [input_format, format_description] = get_file_format(infile);
     const auto tile_id_position = get_tile_info(infile);
     const bool has_tiles = (tile_id_position != 0);
+    const bool has_quals =
+      (input_format == falco::file_format::fastq) || bam_has_quals(infile);
     const auto [n_reads_est, read_len_est, filesize] = [&] {
       // clang-format off
       using falco::file_format;
@@ -231,6 +238,8 @@ get_file_info(const auto &infiles) {
       .size = filesize,
       .n_reads_est = n_reads_est,
       .read_len_est = read_len_est,
+      .has_quals = has_quals,
+      .encoding = falco::encoding::unknown,
       .has_tiles = has_tiles,
       .tile_id_position = tile_id_position,
     });
@@ -309,6 +318,7 @@ main(int argc, char *argv[]) {
     int do_bisulfite{};
     int do_preseq{};
     int do_original_dups{};
+    int do_quals{};
 
     std::uint32_t n_threads{1};
     std::uint32_t max_read_length{};
@@ -437,6 +447,9 @@ main(int argc, char *argv[]) {
     app.add_flag("--tiles,!--no-tiles", do_tiles,
                  "Toggle per-tile quality analysis (default: on)")
       ->option_text(" ");
+    app.add_flag("--quals,!--no-quals", do_quals,
+                 "Toggle quality score analysis (default: on)")
+      ->option_text(" ");
     app.add_flag("--dups,!--no-dups", do_dup_analysis,
                  "Toggle sequence duplication analysis (default: on)")
       ->excludes(orig_dups_opt)
@@ -484,6 +497,7 @@ main(int argc, char *argv[]) {
     mode.set_do_bisulfite(do_bisulfite);
     mode.set_do_preseq(do_preseq);
     mode.set_do_original_dups(do_original_dups);
+    mode.set_do_quals(do_quals);
     mode.set_unassigned();
 
     const auto outdirs = make_outdirs(infiles, outdir, mode.do_stdin());
@@ -538,10 +552,8 @@ main(int argc, char *argv[]) {
                    size_to_units(max_read_length, "bp"));
       std::println("Analyses\n{}", mode.string_verbose());
       std::println("Input files");
-      std::ranges::for_each(infos, [](const auto &info) {
-        std::println("{}\t{}\t{}", info.name, info.description,
-                     size_to_units(info.size, std::string{}));
-      });
+      std::ranges::for_each(
+        infos, [](const auto &info) { std::println("{}", info.to_line()); });
       std::println();
     }
 
diff --git a/src/falco_gc_summary.cpp b/src/falco_gc_summary.cpp
new file mode 100644
index 0000000..e3baf29
--- /dev/null
+++ b/src/falco_gc_summary.cpp
@@ -0,0 +1,160 @@
+// SPDX-License-Identifier: MIT; Copyright 2026 Andrew D Smith
+
+#include "falco_gc_summary.hpp"
+
+#include "falco_utils.hpp"
+
+#include <algorithm>
+#include <array>
+#include <cassert>
+#include <cmath>
+#include <compare>
+#include <cstdint>
+#include <functional>
+#include <iterator>
+#include <numeric>
+#include <ranges>
+#include <tuple>
+#include <utility>
+#include <vector>
+
+[[nodiscard]] auto
+get_theoretical_distribution(const std::vector<double> &gc,
+                             const std::uint64_t total_count)
+  -> std::vector<double> {
+  static constexpr auto mode_width = 0.90;  // ADS: where does this come from?
+  // calculate deviation of a hist from a "normal" with same mode and sd
+  const auto n_bins = std::size(gc);
+  const auto gc_beg = std::cbegin(gc);
+  if (total_count <= 1)  // we will be dividing by (total_count - 1)
+    return std::vector<double>(std::size(gc));
+
+  // get mode
+  const auto mode_itr = std::ranges::max_element(gc);
+  const std::uint64_t mode_pos = std::distance(gc_beg, mode_itr);
+  const auto mode_val = static_cast<double>(*mode_itr);
+
+  // ADS: in case mode is not sharp average nearby values (not clear on why)
+  const auto gt_cut = [&](const double x) { return x < mode_width * mode_val; };
+  const auto right_itr = std::find_if(mode_itr, std::cend(gc), gt_cut);
+  const auto left_itr =
+    std::find_if(std::reverse_iterator(mode_itr), std::crend(gc), gt_cut);
+  const auto n = std::distance(std::reverse_iterator(right_itr), left_itr);
+  const auto mode = static_cast<double>(mode_pos) + as_frac(n - 1, 2.0);
+
+  // theoretical distribution
+  const auto cntr_sq = [m = mode](const auto v) { return (v - m) * (v - m); };
+  const auto sd_term = [&](const auto &x) {
+    return cntr_sq(std::get<0>(x)) * std::get<1>(x);
+  };
+  const auto id_gc =
+    falco::views::enumerate(gc) | std::views::transform(sd_term);
+  const auto sd = std::sqrt(as_frac(
+    std::reduce(std::cbegin(id_gc), std::cend(id_gc)), total_count - 1));
+  const auto to_normal = [&](const auto val) {
+    return std::exp(-cntr_sq(val) / (2.0 * sd * sd));
+  };
+  auto normed = std::views::iota(0U, n_bins) |
+                std::views::transform(to_normal) |
+                std::ranges::to<std::vector>();
+  const auto denom = std::reduce(std::cbegin(normed), std::cend(normed));
+  std::ranges::transform(normed, std::begin(normed),
+                         [&](const auto x) { return x * total_count / denom; });
+  return normed;
+}
+
+[[nodiscard]] auto
+sum_deviation_from_normal(const std::vector<double> &gc) -> double {
+  const auto gc_beg = std::cbegin(gc);
+  const auto gc_end = std::cend(gc);
+  const auto total_count = std::reduce(gc_beg, gc_end);
+  if (total_count <= 1)  // we will be dividing by (total_count - 1)
+    return 0.0;
+  const auto theor =
+    get_theoretical_distribution(gc, static_cast<std::uint64_t>(total_count));
+  const auto diff = [](const auto a, const auto b) { return std::fabs(a - b); };
+  const auto r = std::transform_reduce(gc_beg, gc_end, std::cbegin(theor), 0.0,
+                                       std::plus{}, diff);
+  return pct(as_frac(r, total_count));
+}
+
+[[nodiscard]] auto
+smooth_gc_content(const std::vector<double> &data,
+                  const std::int64_t window_size) -> std::vector<double> {
+  const auto get_mean = [&](std::ranges::viewable_range auto &&r) {
+    return as_frac(std::reduce(std::cbegin(r), std::cend(r)), std::size(r));
+  };
+  assert(std::cmp_less(window_size, std::size(data)));
+  std::vector<double> smoothed;
+  for (auto w = 1; w < (window_size + 1) / 2; ++w)
+    smoothed.push_back(get_mean(
+      std::ranges::subrange(std::cbegin(data), std::cbegin(data) + w)));
+#if __cpp_lib_ranges_slide
+  for (const auto &window : data | std::views::slide(window_size))
+    smoothed.push_back(get_mean(window));
+#else
+  const auto lim = std::cend(data) - window_size;
+  for (auto itr = std::cbegin(data); itr != lim; ++itr)
+    smoothed.push_back(get_mean(std::ranges::subrange(itr, itr + window_size)));
+#endif
+  const auto d_end = std::cend(data);
+  for (auto w = (window_size + 1) / 2; w > 1; --w)
+    smoothed.push_back(get_mean(std::ranges::subrange(d_end - w + 1, d_end)));
+  return smoothed;
+}
+
+[[nodiscard]] static auto
+get_last_bin(const std::int64_t histogram_size, const double next_percent)
+  -> std::int64_t {
+  const auto one_past_last = std::min(
+    histogram_size, static_cast<std::int64_t>(std::ceil(next_percent)));
+  return one_past_last - 1;
+}
+
+static auto
+fill_histogram(
+  const falco::gc_content_t &gc,
+  const std::int32_t gc_idx,
+  const double mass,
+  std::vector<double> &hist  // cppcheck-suppress constParameterReference
+  ) -> void {
+  // ADS: below, not sure best way to do this for all edge cases
+  assert(!hist.empty());
+  const auto curr_percent = gc_idx * mass;
+  const auto next_percent = (gc_idx + 1) * mass;
+  const auto first_bin = static_cast<std::int64_t>(std::floor(curr_percent));
+  const auto last_bin = get_last_bin(std::ssize(hist), next_percent);
+  assert(std::cmp_less(last_bin, std::size(hist)));
+  const auto has_splits = first_bin != last_bin;
+  const auto mass_left =
+    has_splits ? static_cast<double>(first_bin) + 1.0 - curr_percent : mass;
+  const auto mass_right =
+    has_splits ? next_percent - static_cast<double>(last_bin) : mass;
+  const auto hist_begin = std::begin(hist) + first_bin;
+  const auto hist_end = std::cbegin(hist) + last_bin;
+  const auto gc_val = as_frac(gc[gc_idx], mass);
+  for (auto hist_itr = hist_begin; hist_itr <= hist_end; ++hist_itr)
+    *hist_itr += gc_val * ((hist_itr == hist_begin) ? mass_left
+                           : (hist_itr == hist_end) ? mass_right
+                                                    : 1.0);
+}
+
+[[nodiscard]] auto
+combine_gc_content_for_lengths(const std::vector<falco::gc_content_t> &gcs)
+  -> std::vector<double> {
+  static constexpr auto histogram_size = 101;
+  std::vector<double> hist(histogram_size);
+  // ADS: iterate over each gc content vector (per length or length group)
+  for (const auto &gc : gcs) {
+    // ignore any length groups that have no reads
+    if (std::ranges::none_of(gc, [](const auto x) { return x > 0.0; }))
+      continue;
+    // 'mass' is the fraction of each bin in the universal gc-content
+    // histogram that is contributed by each full bin in the current
+    // length-specific gc-content histogram
+    const auto mass = as_frac(histogram_size, std::size(gc));
+    for (auto gc_idx = 0; gc_idx < std::ssize(gc); ++gc_idx)
+      fill_histogram(gc, gc_idx, mass, hist);
+  }
+  return hist;
+}
diff --git a/src/falco_gc_summary.hpp b/src/falco_gc_summary.hpp
new file mode 100644
index 0000000..4a22179
--- /dev/null
+++ b/src/falco_gc_summary.hpp
@@ -0,0 +1,27 @@
+// SPDX-License-Identifier: MIT; Copyright 2026 Andrew D Smith
+
+#ifndef SRC_FALCO_GC_SUMMARY_HPP_
+#define SRC_FALCO_GC_SUMMARY_HPP_
+
+#include "falco_utils.hpp"
+
+#include <cstdint>
+#include <vector>
+
+[[nodiscard]] auto
+combine_gc_content_for_lengths(const std::vector<falco::gc_content_t> &gcs)
+  -> std::vector<double>;
+
+[[nodiscard]] auto
+smooth_gc_content(const std::vector<double> &data,
+                  const std::int64_t window_size) -> std::vector<double>;
+
+[[nodiscard]] auto
+get_theoretical_distribution(const std::vector<double> &gc,
+                             const std::uint64_t total_count)
+  -> std::vector<double>;
+
+[[nodiscard]] auto
+sum_deviation_from_normal(const std::vector<double> &gc) -> double;
+
+#endif  // SRC_FALCO_GC_SUMMARY_HPP_
diff --git a/src/falco_grade.cpp b/src/falco_grade.cpp
index 61295fa..28896e0 100644
--- a/src/falco_grade.cpp
+++ b/src/falco_grade.cpp
@@ -1,6 +1,7 @@
 // SPDX-License-Identifier: MIT; Copyright 2026 Andrew D Smith
 
 #include "falco_grade.hpp"
+#include "falco_gc_summary.hpp"
 
 #include <algorithm>
 #include <cassert>
diff --git a/src/falco_utils.cpp b/src/falco_utils.cpp
index da425d7..8ab65fb 100644
--- a/src/falco_utils.cpp
+++ b/src/falco_utils.cpp
@@ -2,19 +2,13 @@
 
 #include "falco_utils.hpp"
 
-#include <algorithm>
-#include <array>
 #include <cmath>
-#include <compare>
 #include <cstdint>
-#include <ctime>  // for std::localtime
+#include <ctime>
 #include <format>
-#include <iomanip>  // for std::put_time
+#include <iomanip>
 #include <sstream>
 #include <string>
-#include <tuple>
-#include <utility>
-#include <vector>
 
 [[nodiscard]] auto
 size_to_units(const std::int64_t s, const std::string &suffix) -> std::string {
@@ -30,129 +24,6 @@ size_to_units(const std::int64_t s, const std::string &suffix) -> std::string {
   return std::format("{}", s);
 };
 
-[[nodiscard]] auto
-get_theoretical_distribution(const std::vector<double> &gc,
-                             const std::uint64_t total_count)
-  -> std::vector<double> {
-  static constexpr auto mode_width = 0.90;  // ADS: where does this come from?
-  // calculate deviation of a hist from a "normal" with same mode and sd
-  const auto n_bins = std::size(gc);
-  const auto gc_beg = std::cbegin(gc);
-  if (total_count <= 1)  // we will be dividing by (total_count - 1)
-    return std::vector<double>(std::size(gc));
-
-  // get mode
-  const auto mode_itr = std::ranges::max_element(gc);
-  const std::uint64_t mode_pos = std::distance(gc_beg, mode_itr);
-  const auto mode_val = static_cast<double>(*mode_itr);
-
-  // ADS: in case mode is not sharp average nearby values (not clear on why)
-  const auto gt_cut = [&](const double x) { return x < mode_width * mode_val; };
-  const auto right_itr = std::find_if(mode_itr, std::cend(gc), gt_cut);
-  const auto left_itr =
-    std::find_if(std::reverse_iterator(mode_itr), std::crend(gc), gt_cut);
-  const auto n = std::distance(std::reverse_iterator(right_itr), left_itr);
-  const auto mode = static_cast<double>(mode_pos) + as_frac(n - 1, 2.0);
-
-  // theoretical distribution
-  const auto cntr_sq = [m = mode](const auto v) { return (v - m) * (v - m); };
-  const auto sd_term = [&](const auto &x) {
-    return cntr_sq(std::get<0>(x)) * std::get<1>(x);
-  };
-  const auto id_gc =
-    falco::views::enumerate(gc) | std::views::transform(sd_term);
-  const auto sd = std::sqrt(as_frac(
-    std::reduce(std::cbegin(id_gc), std::cend(id_gc)), total_count - 1));
-  const auto to_normal = [&](const auto val) {
-    return std::exp(-cntr_sq(val) / (2.0 * sd * sd));
-  };
-  auto normed = std::views::iota(0U, n_bins) |
-                std::views::transform(to_normal) |
-                std::ranges::to<std::vector>();
-  const auto denom = std::reduce(std::cbegin(normed), std::cend(normed));
-  std::ranges::transform(normed, std::begin(normed),
-                         [&](const auto x) { return x * total_count / denom; });
-  return normed;
-}
-
-[[nodiscard]] auto
-sum_deviation_from_normal(const std::vector<double> &gc) -> double {
-  const auto gc_beg = std::cbegin(gc);
-  const auto gc_end = std::cend(gc);
-  const auto total_count = std::reduce(gc_beg, gc_end);
-  if (total_count <= 1)  // we will be dividing by (total_count - 1)
-    return 0.0;
-  const auto theor =
-    get_theoretical_distribution(gc, static_cast<std::uint64_t>(total_count));
-  const auto diff = [](const auto a, const auto b) { return std::fabs(a - b); };
-  const auto r = std::transform_reduce(gc_beg, gc_end, std::cbegin(theor), 0.0,
-                                       std::plus{}, diff);
-  return pct(as_frac(r, total_count));
-}
-
-[[nodiscard]] auto
-smooth_gc_content(const std::vector<double> &data,
-                  const std::int64_t window_size) -> std::vector<double> {
-  const auto get_mean = [&](std::ranges::viewable_range auto &&r) {
-    return as_frac(std::reduce(std::cbegin(r), std::cend(r)), std::size(r));
-  };
-  assert(std::cmp_less(window_size, std::size(data)));
-  std::vector<double> smoothed;
-  for (auto w = 1; w < (window_size + 1) / 2; ++w)
-    smoothed.push_back(get_mean(
-      std::ranges::subrange(std::cbegin(data), std::cbegin(data) + w)));
-#if __cpp_lib_ranges_slide
-  for (const auto &window : data | std::views::slide(window_size))
-    smoothed.push_back(get_mean(window));
-#else
-  const auto lim = std::cend(data) - window_size;
-  for (auto itr = std::cbegin(data); itr != lim; ++itr)
-    smoothed.push_back(get_mean(std::ranges::subrange(itr, itr + window_size)));
-#endif
-  const auto d_end = std::cend(data);
-  for (auto w = (window_size + 1) / 2; w > 1; --w)
-    smoothed.push_back(get_mean(std::ranges::subrange(d_end - w + 1, d_end)));
-  return smoothed;
-}
-
-[[nodiscard]] auto
-combine_gc_content_for_lengths(const std::vector<falco::gc_content_t> &gcs)
-  -> std::vector<double> {
-  static constexpr auto histogram_size = 101;
-  std::vector<double> hist(histogram_size);
-  for (const auto &gc : gcs) {
-    if (std::ranges::none_of(gc, [](const auto x) { return x > 0.0; }))
-      continue;
-    const auto increment = as_frac(histogram_size, std::size(gc));
-    for (auto gc_idx = 0U; gc_idx < std::size(gc); ++gc_idx) {
-      const auto curr_percent = gc_idx * increment;
-      const auto next_percent = (gc_idx + 1) * increment;
-      const auto start_in_hist =
-        static_cast<std::int64_t>(std::floor(curr_percent));
-      // ADS: below, not sure best way to do this for all edge cases
-      const auto stop_in_hist =
-        static_cast<std::int64_t>(std::min(static_cast<double>(histogram_size),
-                                           std::ceil(next_percent))) -
-        1;
-      assert(stop_in_hist < histogram_size);
-      const auto splits = start_in_hist != stop_in_hist;
-      const auto frac_left =
-        splits ? static_cast<double>(start_in_hist) + 1.0 - curr_percent
-               : increment;
-      const auto frac_right =
-        splits ? next_percent - static_cast<double>(stop_in_hist) : increment;
-      const auto hist_begin = std::begin(hist) + start_in_hist;
-      const auto hist_end = std::cbegin(hist) + stop_in_hist;
-      const auto gc_val = as_frac(gc[gc_idx], increment);
-      for (auto hist_itr = hist_begin; hist_itr <= hist_end; ++hist_itr)
-        *hist_itr += gc_val * ((hist_itr == hist_begin) ? frac_left
-                               : (hist_itr == hist_end) ? frac_right
-                                                        : 1.0);
-    }
-  }
-  return hist;
-}
-
 [[nodiscard]] auto
 get_program_start_time() -> std::chrono::time_point<std::chrono::system_clock> {
   static const auto start_time = std::chrono::system_clock::now();
diff --git a/src/falco_utils.hpp b/src/falco_utils.hpp
index 504fa4d..aaf934d 100644
--- a/src/falco_utils.hpp
+++ b/src/falco_utils.hpp
@@ -106,22 +106,6 @@ resize_gc_content(const std::uint32_t updated_length,
     gc.resize(prev_size + idx + 1);
 }
 
-[[nodiscard]] auto
-combine_gc_content_for_lengths(const std::vector<falco::gc_content_t> &gcs)
-  -> std::vector<double>;
-
-[[nodiscard]] auto
-smooth_gc_content(const std::vector<double> &data,
-                  const std::int64_t window_size) -> std::vector<double>;
-
-[[nodiscard]] auto
-get_theoretical_distribution(const std::vector<double> &gc,
-                             const std::uint64_t total_count)
-  -> std::vector<double>;
-
-[[nodiscard]] auto
-sum_deviation_from_normal(const std::vector<double> &gc) -> double;
-
 [[nodiscard]] auto
 get_run_duration(const auto start_time) {
   using namespace std::literals::chrono_literals;
diff --git a/src/fastq_bgzf_file.cpp b/src/fastq_bgzf_file.cpp
index acd5b6d..edf6de2 100644
--- a/src/fastq_bgzf_file.cpp
+++ b/src/fastq_bgzf_file.cpp
@@ -99,23 +99,23 @@ fastq_bgzf_file::load_next(const std::int32_t file_id,
     input_last = n_to_keep;
   }
   br.release();  // use like monotonic_buffer_resource
-  // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic)
-  auto in_itr = std::data(input_buffer) + input_last;
+  auto in_itr = std::begin(input_buffer) + input_last;
   while (br.task_ready() && has_in()) {
-    auto task = br.get_decomp_task(in_itr);
+    auto task = br.get_decomp_task(std::to_address(in_itr));
     // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic)
     in_itr += task.size;
     ++n_tasks;
     tq.push(file_id, std::move(task));
   }
-  input_last = std::distance(std::data(input_buffer), in_itr);
+  input_last = std::distance(std::begin(input_buffer), in_itr);
 }
 
 auto
-fastq_bgzf_file::get_chunks(const std::int64_t n_chunks,
+fastq_bgzf_file::get_chunks(std::int64_t n_chunks,
                             const std::int32_t file_id,
                             task_queue &tq,
                             std::atomic_int32_t &n_tasks) -> void {
+  static constexpr std::int64_t min_chunk_size = 64 * 1024;
   static constexpr auto rec_lines = 4;  // FASTQ
   assert(n_chunks > 0);
   std::swap(input_buffer, output_buffer);
@@ -142,6 +142,8 @@ fastq_bgzf_file::get_chunks(const std::int64_t n_chunks,
   };
   // NOLINTEND(cppcoreguidelines-pro-bounds-pointer-arithmetic)
   const auto n_bytes_available = output_last;
+  n_chunks = std::min((n_bytes_available + min_chunk_size - 1) / min_chunk_size,
+                      n_chunks);
   const auto [chunk_size, remainder] = std::div(n_bytes_available, n_chunks);
   assert(n_chunks > 0);
   std::int64_t start_pos{};
diff --git a/src/fastq_bgzf_file.hpp b/src/fastq_bgzf_file.hpp
index 064274d..0023291 100644
--- a/src/fastq_bgzf_file.hpp
+++ b/src/fastq_bgzf_file.hpp
@@ -18,7 +18,7 @@
 struct task_queue;
 
 class fastq_bgzf_file {
-  static constexpr std::int64_t min_buf_size = 64 * 1024;
+  static constexpr std::int64_t min_buf_size = 256 * 1024;
   std::vector<char> input_buffer;
   std::vector<char> output_buffer;
   std::int64_t input_last{};
@@ -86,17 +86,17 @@ class fastq_bgzf_file {
     br.release();
   }
 
-  auto
-  get_chunks(const std::int64_t n_chunks,
-             const std::int32_t file_id,
-             task_queue &tq,
-             std::atomic_int32_t &n_tasks) -> void;
-
   auto
   load_next(const std::int32_t file_id,
             task_queue &tq,
             std::atomic_int32_t &n_tasks) -> void;
 
+  auto
+  get_chunks(std::int64_t n_chunks,
+             const std::int32_t file_id,
+             task_queue &tq,
+             std::atomic_int32_t &n_tasks) -> void;
+
   auto
   make_tasks_inflate(const std::int32_t file_id,
                      task_queue &tq,
diff --git a/src/fastq_file.cpp b/src/fastq_file.cpp
index b564851..2abd3e7 100644
--- a/src/fastq_file.cpp
+++ b/src/fastq_file.cpp
@@ -100,7 +100,9 @@ fastq_file::reset() -> void {
 }
 
 static inline auto
-mmap_fastq(const int fd, const std::int64_t offset, const std::int64_t length,
+mmap_fastq(const int fd,
+           const std::int64_t offset,
+           const std::int64_t length,
            auto &data) {
   static constexpr auto prot = PROT_READ;
   static constexpr auto flags = MAP_PRIVATE;
@@ -130,10 +132,11 @@ fastq_file::load_next() -> void {
 }
 
 auto
-fastq_file::get_chunks(const std::int64_t n_chunks,
+fastq_file::get_chunks(std::int64_t n_chunks,
                        const std::int32_t file_id,
                        task_queue &tq,
                        std::atomic_int32_t &n_tasks) -> void {
+  static constexpr std::int64_t min_chunk_size = 64 * 1024;
   static constexpr auto rec_lines = 4;  // FASTQ
   assert(n_chunks > 0);
   const auto beg_itr = std::begin(buffer);
@@ -154,7 +157,9 @@ fastq_file::get_chunks(const std::int64_t n_chunks,
     return pos;
   };
   // clang-format on
-  const auto [chunk_size, remainder] = std::ldiv(std::ssize(buffer), n_chunks);
+  n_chunks = std::min(
+    (std::ssize(buffer) + min_chunk_size - 1) / min_chunk_size, n_chunks);
+  auto [chunk_size, remainder] = std::ldiv(std::ssize(buffer), n_chunks);
   auto start_pos = beg_itr;
   auto chunk_end = start_pos;
   for (const auto chunk_idx : std::views::iota(0, n_chunks)) {
@@ -171,3 +176,14 @@ fastq_file::get_chunks(const std::int64_t n_chunks,
   }
   last = chunk_end;
 }
+
+auto
+fastq_file::make_tasks(const std::int64_t n_threads,
+                       const std::int32_t file_id,
+                       task_queue &tq,
+                       std::atomic_int32_t &n_tasks) -> void {
+  static constexpr auto n_chunks_per_thread = 8;
+  const auto n_chunks = n_chunks_per_thread * n_threads;
+  load_next();
+  get_chunks(n_chunks, file_id, tq, n_tasks);
+}
diff --git a/src/fastq_file.hpp b/src/fastq_file.hpp
index 2d0763c..19d5ec1 100644
--- a/src/fastq_file.hpp
+++ b/src/fastq_file.hpp
@@ -15,7 +15,7 @@
 struct task_queue;
 
 struct fastq_file {
-  static constexpr auto min_buf_size = 65536;
+  static constexpr auto min_buf_size = 256 * 1024;
   std::int64_t target_length{};
   std::int64_t filesize{};
   char *mmap_data{};
@@ -46,9 +46,6 @@ struct fastq_file {
     last{src.last}                     //
   {}
 
-  auto
-  reset() -> void;
-
   ~fastq_file() {
     reset();
     close(fd);  // will always have been opened using a filename
@@ -56,11 +53,31 @@ struct fastq_file {
 
   operator bool() const { return stop_in_file != filesize; }
 
+  friend auto
+  reset(fastq_file &reads_file) -> void;
+
+  friend auto
+  make_tasks(fastq_file &reads_file,
+             const std::int64_t n_threads,
+             const std::int32_t file_id,
+             task_queue &tq,
+             std::atomic_int32_t &n_tasks) -> void;
+
+private:
+  auto
+  reset() -> void;
+
   auto
   load_next() -> void;
 
   auto
-  get_chunks(const std::int64_t n_chunks,
+  get_chunks(std::int64_t n_chunks,
+             const std::int32_t file_id,
+             task_queue &tq,
+             std::atomic_int32_t &n_tasks) -> void;
+
+  auto
+  make_tasks(const std::int64_t n_threads,
              const std::int32_t file_id,
              task_queue &tq,
              std::atomic_int32_t &n_tasks) -> void;
@@ -76,11 +93,8 @@ make_tasks(fastq_file &reads_file,
            const std::int32_t file_id,
            task_queue &tq,
            std::atomic_int32_t &n_tasks) -> void {
-  static constexpr auto n_chunks_per_thread = 8;
-  const auto n_chunks = n_chunks_per_thread * n_threads;
   n_tasks = 1;  // for current task, which makes tasks
-  reads_file.load_next();
-  reads_file.get_chunks(n_chunks, file_id, tq, n_tasks);
+  reads_file.make_tasks(n_threads, file_id, tq, n_tasks);
 }
 
 inline auto
diff --git a/src/fastq_gz_file.cpp b/src/fastq_gz_file.cpp
index 67fdf06..e87e104 100644
--- a/src/fastq_gz_file.cpp
+++ b/src/fastq_gz_file.cpp
@@ -50,11 +50,11 @@ estimate_n_reads_fastq_gz(const std::string &filename)
 }
 
 auto
-fastq_gz_file::get_chunks(const std::int64_t n_chunks,  //
-                          const std::int32_t file_id,   //
-                          task_queue &tq,               //
-                          std::atomic_int32_t &n_tasks  //
-                          ) -> void {
+fastq_gz_file::get_chunks(std::int64_t n_chunks,
+                          const std::int32_t file_id,
+                          task_queue &tq,
+                          std::atomic_int32_t &n_tasks) -> void {
+  static constexpr std::int64_t min_chunk_size = 64 * 1024;
   static constexpr auto rec_lines = 4;  // FASTQ
   assert(n_chunks > 0);
   const auto data = std::data(outbuf);
@@ -77,8 +77,9 @@ fastq_gz_file::get_chunks(const std::int64_t n_chunks,  //
     return p;
   };
   // NOLINTEND(cppcoreguidelines-pro-bounds-pointer-arithmetic)
-  const auto n_bytes_available =
-    buf_sz - cursor;  // ADS: I think cursor must always be 0 here
+  const auto n_bytes_available = buf_sz - cursor;  // ADS: cursor always 0 here?
+  n_chunks = std::min((n_bytes_available + min_chunk_size - 1) / min_chunk_size,
+                      n_chunks);
   const auto [chunk_size, remainder] = std::div(n_bytes_available, n_chunks);
   assert(n_chunks > 0);
   std::int64_t start_pos = cursor;  // ADS: I think cursor must always be 0 here
diff --git a/src/fastq_gz_file.hpp b/src/fastq_gz_file.hpp
index 40bf782..f4f977e 100644
--- a/src/fastq_gz_file.hpp
+++ b/src/fastq_gz_file.hpp
@@ -37,7 +37,7 @@ struct task_queue;
 #ifdef HAVE_ISAL
 
 class fastq_gz_file {
-  static constexpr auto min_buf_size = 64 * 1024;
+  static constexpr auto min_buf_size = 256 * 1024;
   static constexpr auto inflate_err_msg =
     R"(Failure during decompression by ISAL. Error code is {}.  Please check that the
 input file is not corrupted by decompressing with gunzip. If the input file is
@@ -177,7 +177,7 @@ not corrupted please file a bug report at the Falco repo on GitHub.)";
 #else  // use bgzf for ordinary gz files
 
 class fastq_gz_file {
-  static constexpr auto min_buf_size = 64 * 1024;
+  static constexpr auto min_buf_size = 256 * 1024;
   char *buf_data{};         // start of the outbuf
   std::int64_t buf_sz{};    // amount of output buffer used
   std::int64_t buf_size{};  // size of allocated buffer
diff --git a/src/fastq_stdin.cpp b/src/fastq_stdin.cpp
index 39f2af0..df935cb 100644
--- a/src/fastq_stdin.cpp
+++ b/src/fastq_stdin.cpp
@@ -13,6 +13,7 @@
 #include <cstdlib>
 #include <memory>
 #include <ranges>
+#include <span>
 #include <stdexcept>
 #include <string>
 #include <system_error>
@@ -33,10 +34,11 @@ estimate_n_reads_fastq_stdin(const std::string &)
 }
 
 auto
-fastq_stdin::get_chunks(const std::int64_t n_chunks,
+fastq_stdin::get_chunks(std::int64_t n_chunks,
                         const std::int32_t file_id,
                         task_queue &tq,
                         std::atomic_int32_t &n_tasks) -> void {
+  static constexpr std::int64_t min_chunk_size = 64 * 1024;
   static constexpr auto rec_lines = 4;  // FASTQ
   const auto beg_itr = std::begin(buffer);
   const auto end_itr = last;
@@ -58,6 +60,8 @@ fastq_stdin::get_chunks(const std::int64_t n_chunks,
   };
   // clang-format on
   const std::int64_t n_bytes_available = std::distance(beg_itr, last);
+  n_chunks = std::min((n_bytes_available + min_chunk_size - 1) / min_chunk_size,
+                      n_chunks);
   const auto [chunk_size, remainder] = std::div(n_bytes_available, n_chunks);
   auto start_itr = beg_itr;
   auto chunk_end = start_itr;
@@ -105,7 +109,7 @@ validate_fastq(const auto &buffer) -> bool {
     prev_nl = (c == '\n');
     n_lines += prev_nl;
     if (++n_bytes == n_bytes_to_validate)
-      return true;
+      break;
   }
   return true;
 }
@@ -125,14 +129,15 @@ fastq_stdin::load_next() -> void {
 }
 
 auto
-fastq_stdin::make_tasks(const std::int64_t n_chunks,  //
-                        const std::int32_t file_id,   //
-                        task_queue &tq,               //
+fastq_stdin::make_tasks(const std::int64_t n_threads,
+                        const std::int32_t file_id,
+                        task_queue &tq,
                         std::atomic_int32_t &n_tasks) -> void {
-  n_tasks = 1;  // for current task, which makes more tasks
+  static constexpr auto n_chunks_per_thread = 8;
+  const auto n_chunks = n_chunks_per_thread * n_threads;
   shift_output_buffer();
   load_next();
-  if (!validate_fastq(buffer))
+  if (!validate_fastq(std::span(std::cbegin(buffer), last)))
     throw std::runtime_error("input appears not to be FASTQ");
   get_chunks(n_chunks, file_id, tq, n_tasks);
 }
diff --git a/src/fastq_stdin.hpp b/src/fastq_stdin.hpp
index 1603b08..c6060d8 100644
--- a/src/fastq_stdin.hpp
+++ b/src/fastq_stdin.hpp
@@ -14,7 +14,7 @@
 struct task_queue;
 
 struct fastq_stdin {
-  static constexpr auto min_buf_size = 64 * 1024;
+  static constexpr auto min_buf_size = 256 * 1024;
   std::vector<char> buffer;
   std::vector<char>::iterator cursor;
   std::vector<char>::iterator last;
@@ -31,6 +31,17 @@ struct fastq_stdin {
   ~fastq_stdin() = default;
   // clang-format on
 
+  friend auto
+  reset(fastq_stdin &reads_file) -> void;
+
+  friend auto
+  make_tasks(fastq_stdin &reads_file,
+             const std::int64_t n_threads,
+             const std::int32_t file_id,
+             task_queue &tq,
+             std::atomic_int32_t &n_tasks) -> void;
+
+private:
   auto
   reset() -> void {
     buffer.clear();
@@ -40,14 +51,13 @@ struct fastq_stdin {
   }
 
   auto
-  make_tasks(const std::int64_t n_chunks,
+  make_tasks(const std::int64_t n_threads,
              const std::int32_t file_id,
              task_queue &tq,
              std::atomic_int32_t &n_tasks) -> void;
 
-private:
   auto
-  get_chunks(const std::int64_t n_chunks,
+  get_chunks(std::int64_t n_chunks,
              const std::int32_t file_id,
              task_queue &tq,
              std::atomic_int32_t &n_tasks) -> void;
@@ -64,14 +74,13 @@ estimate_n_reads_fastq_stdin(const std::string &filename)
   -> std::tuple<std::uint64_t, std::uint64_t, std::int64_t>;
 
 inline auto
-make_tasks(fastq_stdin &reads_file,       //
-           const std::int64_t n_threads,  //
-           const std::int32_t file_id,    //
-           task_queue &tq,                //
+make_tasks(fastq_stdin &reads_file,
+           const std::int64_t n_threads,
+           const std::int32_t file_id,
+           task_queue &tq,
            std::atomic_int32_t &n_tasks) -> void {
-  static constexpr auto n_chunks_per_thread = 8;
-  const auto n_chunks = n_chunks_per_thread * n_threads;
-  reads_file.make_tasks(n_chunks, file_id, tq, n_tasks);
+  n_tasks = 1;  // for current task, which makes more tasks
+  reads_file.make_tasks(n_threads, file_id, tq, n_tasks);
 }
 
 inline auto
diff --git a/src/file_info.cpp b/src/file_info.cpp
new file mode 100644
index 0000000..6c903df
--- /dev/null
+++ b/src/file_info.cpp
@@ -0,0 +1,15 @@
+// SPDX-License-Identifier: MIT; Copyright 2026 Andrew D Smith
+
+#include "file_info.hpp"
+#include "falco_utils.hpp"
+
+#include <format>
+#include <string>
+
+[[nodiscard]] auto
+file_info::to_line() const -> std::string {
+  return std::format("{}\t{}\t{}", name,
+                     description + (!has_quals ? std::string(" (quals missing)")
+                                               : std::string{}),
+                     size_to_units(size, std::string{}));
+}
diff --git a/src/file_info.hpp b/src/file_info.hpp
index c8c6415..d99b23f 100644
--- a/src/file_info.hpp
+++ b/src/file_info.hpp
@@ -3,14 +3,16 @@
 #ifndef SRC_FILE_INFO_HPP_
 #define SRC_FILE_INFO_HPP_
 
-#include "falco_file_format.hpp"
-#include "quality_score.hpp"
-
 #include "nlohmann/json.hpp"
 
 #include <cstdint>
 #include <string>
 
+namespace falco {
+enum class encoding : std::uint8_t;
+enum class file_format : std::uint8_t;
+}  // namespace falco
+
 struct file_info {
   std::string name;
   falco::file_format format{};
@@ -18,16 +20,20 @@ struct file_info {
   std::int64_t size{};
   std::uint64_t n_reads_est{};
   std::uint64_t read_len_est{};
+  bool has_quals{};
   falco::encoding encoding{};
   bool has_tiles{};
   std::uint32_t tile_id_position{};
 
   [[nodiscard]] auto
-  string() const -> std::string {
+  to_string() const -> std::string {
     static constexpr auto n_indent = 4;
     return nlohmann::json(*this).dump(n_indent);
   }
 
+  [[nodiscard]] auto
+  to_line() const -> std::string;
+
   auto
   set_encoding(const falco::encoding &e) {
     encoding = e;
@@ -40,6 +46,7 @@ struct file_info {
                                  size,
                                  n_reads_est,
                                  read_len_est,
+                                 has_quals,
                                  encoding,
                                  has_tiles,
                                  tile_id_position);
diff --git a/src/fqrec.hpp b/src/fqrec.hpp
index 06185a1..6737415 100644
--- a/src/fqrec.hpp
+++ b/src/fqrec.hpp
@@ -36,7 +36,7 @@ struct fqrec {
   size() const { return static_cast<std::int32_t>(std::distance(r, o)) - 1; }
   [[nodiscard]] operator bool() const { return n != nullptr; }
   [[nodiscard]] auto
-  string() const -> std::string { return {n, e}; }
+  to_string() const -> std::string { return {n, e}; }
 };
 // NOLINTBEGIN(cppcoreguidelines-pro-bounds-pointer-arithmetic)
 [[nodiscard]] constexpr auto get_name(const fqrec &rec) { return rec.n; }
diff --git a/src/html.cpp b/src/html.cpp
index 7ffd728..cdf5833 100644
--- a/src/html.cpp
+++ b/src/html.cpp
@@ -3,6 +3,7 @@
 #include "html.hpp"
 
 #include "base_groups.hpp"
+#include "falco_gc_summary.hpp"
 #include "falco_grade.hpp"
 #include "falco_utils.hpp"
 #include "file_info.hpp"
@@ -417,7 +418,7 @@ yaxis: {{title: "tile", type: "category"}},
   };
   const auto z = std::views::transform(centered | std::views::values, format1);
 
-  assert(std::size(std::cbegin(centered)->second) == std::size(groups));
+  assert(std::size(std::cbegin(centered)->second) <= std::size(groups));
   const auto grade = grades.grade(label);
   const auto title = grades.get_title(label);
   return fmt::format(
diff --git a/src/quality_score.cpp b/src/quality_score.cpp
index c652205..16ff895 100644
--- a/src/quality_score.cpp
+++ b/src/quality_score.cpp
@@ -2,16 +2,22 @@
 
 #include "quality_score.hpp"
 
+#include "bamrec.hpp"
 #include "falco_file_format.hpp"
 #include "file_info.hpp"
 
+#include <htslib/sam.h>
+
 #include <algorithm>
 #include <cassert>  // IWYU pragma: keep
+#include <cerrno>
 #include <format>
 #include <iterator>
+#include <memory>
 #include <ranges>
 #include <stdexcept>
 #include <string>
+#include <system_error>
 #include <utility>
 #include <vector>
 
@@ -24,9 +30,11 @@ identify_encoding_impl(const auto &qual_counts,
   };
   if (qual_counts.empty())
     return falco::encoding::unknown;
-  const auto min_qual =
-    std::ranges::min(qual_counts | std::views::transform(first_non_zero)) +
-    qual_offset;
+  const auto min_qual_observed =
+    std::ranges::min(qual_counts | std::views::transform(first_non_zero));
+  if (min_qual_observed == std::ssize(qual_counts.front()))
+    return falco::encoding::unknown;
+  const auto min_qual = min_qual_observed + qual_offset;
   if (min_qual > falco::max_qual_val)
     throw std::runtime_error("invalid qual score: " + std::to_string(min_qual));
   if (min_qual < falco::sanger_min_qual)
@@ -73,3 +81,29 @@ encoding_to_string(const falco::encoding e) -> std::string {
   assert(u < std::size(falco::format_labels));
   return falco::format_labels[u];
 }
+
+[[nodiscard]] auto
+bam_has_quals(const std::string &filename) -> bool {
+  static constexpr auto n_records = 100;
+  const std::unique_ptr<samFile, int (*)(samFile *)> in(
+    hts_open(std::data(filename), "r"), &hts_close);
+  if (!in)
+    throw std::runtime_error("failed to open BAM/SAM file: " + filename);
+  const std::unique_ptr<sam_hdr_t, void (*)(sam_hdr_t *)> h(
+    sam_hdr_read(in.get()), &sam_hdr_destroy);
+  if (!h)
+    throw std::system_error(std::make_error_code(std::errc(errno)),
+                            "failed to read header: " + filename);
+  for (auto i = 0; i < n_records; ++i) {
+    std::unique_ptr<bam1_t, void (*)(bam1_t *)> b(bam_init1(), &bam_destroy1);
+    const auto r = sam_read1(in.get(), h.get(), b.get());  // -1 on EOF
+    if (r == -1)  // EOF (htslib/bgzf.h)
+      break;
+    if (r < -1)  // ERROR (htslib/bgzf.h)
+      throw std::runtime_error("failed reading BAM/SAM record: " + filename);
+    if (static_cast<char>(*bam_get_qual(b)) ==
+        bamrec::get_qual_missing_code())  // (htslib/sam.h)
+      return false;
+  }
+  return true;
+}
diff --git a/src/quality_score.hpp b/src/quality_score.hpp
index 205f9b1..831b228 100644
--- a/src/quality_score.hpp
+++ b/src/quality_score.hpp
@@ -72,6 +72,9 @@ adjust_fastq_qual_encoding(std::vector<falco::qual_array> &qual_by_pos,
                            falco::qual_array &qual_by_read,
                            const falco::encoding enc) -> void;
 
+[[nodiscard]] auto
+bam_has_quals(const std::string &filename) -> bool;
+
 [[nodiscard]] static inline auto
 count_quals(auto qual_itr,
             const auto qual_end,
diff --git a/src/results_collector.hpp b/src/results_collector.hpp
index 552b77d..5e521fc 100644
--- a/src/results_collector.hpp
+++ b/src/results_collector.hpp
@@ -71,8 +71,8 @@ struct alignas(assumed_page_size) results_collector {
     }
     else
       dr.initialize(mode, info);
-    do_tiles = mode.do_tiles() && info.has_tiles;
     do_kmers = mode.do_kmers();
+    do_tiles = mode.do_quals() && mode.do_tiles() && info.has_tiles;
     if (do_tiles)
       tp.init(info);
   }
diff --git a/src/results_summary.cpp b/src/results_summary.cpp
index 6699da3..8743fdd 100644
--- a/src/results_summary.cpp
+++ b/src/results_summary.cpp
@@ -3,6 +3,7 @@
 #include "results_summary.hpp"
 
 #include "base_groups.hpp"
+#include "falco_gc_summary.hpp"
 #include "falco_utils.hpp"
 #include "html.hpp"
 #include "report.hpp"
@@ -47,7 +48,7 @@ results_summary::initialize() -> void {
     return a + nuc[guanine_index] + nuc[cytosine_index];
   };
   total_gc = std::accumulate(std::cbegin(base_counts), std::cend(base_counts),
-                             0ul, gc_acc);
+                             0LU, gc_acc);
   median_read_len =
     min_read_len_itr == std::cend(lengths) ? 0LU : median_tabular(lengths);
 
@@ -94,10 +95,10 @@ results_summary::assign_grades() -> void {
     grades.emplace("overrepresented",
                    get_grade_overrepresented(dup_summary.n_reads, dr));
 
-  if (mode.do_qual_base())
+  if (mode.do_qual_base() && info.has_quals)
     grades.emplace("quality_base", get_grade_quality_base(qual_by_pos));
 
-  if (mode.do_qual_seq())
+  if (mode.do_qual_seq() && info.has_quals)
     grades.emplace("quality_sequence",
                    get_grade_quality_sequence(qual_by_read));
 
@@ -148,11 +149,11 @@ results_summary::get_report() const -> std::string {
     sections.emplace("overrepresented",
                      overrepresented_report(overrep, grades));
 
-  if (mode.do_qual_base())
+  if (mode.do_qual_base() && info.has_quals)
     sections.emplace("quality_base",
                      quality_base_report(qual_by_pos, groups, grades));
 
-  if (mode.do_qual_seq())
+  if (mode.do_qual_seq() && info.has_quals)
     sections.emplace("quality_sequence",
                      quality_sequence_report(qual_by_read, grades));
 
@@ -213,11 +214,11 @@ results_summary::get_html() -> std::string {
   if (mode.do_overrep())
     sections.emplace("overrepresented", overrepresented_html(overrep, grades));
 
-  if (mode.do_qual_base())
+  if (mode.do_qual_base() && info.has_quals)
     sections.emplace("quality_base",
                      quality_base_html(qual_by_pos, groups, grades));
 
-  if (mode.do_qual_seq())
+  if (mode.do_qual_seq() && info.has_quals)
     sections.emplace("quality_sequence",
                      quality_sequence_html(qual_by_read, grades));
 
diff --git a/src/run_mode.cpp b/src/run_mode.cpp
index 1713ad8..499967d 100644
--- a/src/run_mode.cpp
+++ b/src/run_mode.cpp
@@ -43,6 +43,7 @@
 // do_bisulfite (changes grading and nothing else)
 // do_original_dups (use the original duplication mode from FastQC and Falco v1)
 // do_preseq (make another output file for input to preseq)
+// do_quals (controls any quality score analysis in case no qual scores present)
 //
 // ADS (2026-09-08 edit): added do_stdin because otherwise it's a floating
 // almost global variable.
@@ -97,6 +98,7 @@ run_mode::set_unassigned() -> void {
   if (do_bisulfite_ == 0) do_bisulfite_ = do_bisulfite_default;
   if (do_preseq_ == 0) do_preseq_ = do_preseq_default;
   if (do_original_dups_ == 0) do_original_dups_ = do_original_dups_default;
+  if (do_quals_ == 0) do_quals_ = do_quals_default;
   //
   if (do_adap_ == 0) do_adap_ = do_adap_default;
   if (do_dups_ == 0) do_dups_ = do_dups_default;
@@ -110,6 +112,11 @@ run_mode::set_unassigned() -> void {
   if (do_length_ == 0) do_length_ = do_length_default;
   if (do_tiles_ == 0) do_tiles_ = do_tiles_default;
   // clang-format on
+
+  if (do_quals_ == -1) {
+    do_qual_seq_ = -1;
+    do_qual_base_ = -1;
+  }
 }
 
 [[nodiscard]] auto
diff --git a/src/run_mode.hpp b/src/run_mode.hpp
index 4e838a1..146a8ce 100644
--- a/src/run_mode.hpp
+++ b/src/run_mode.hpp
@@ -22,6 +22,7 @@ class run_mode {
   [[nodiscard]] auto do_bisulfite() const -> bool { return do_bisulfite_ == 1; }
   [[nodiscard]] auto do_preseq() const -> bool { return do_preseq_ == 1; }
   [[nodiscard]] auto do_original_dups() const -> bool { return do_original_dups_ == 1; }
+  [[nodiscard]] auto do_quals() const -> bool { return do_quals_ == 1; }
   //
   [[nodiscard]] auto do_adap() const -> bool { return do_adap_ == 1; }
   [[nodiscard]] auto do_dups() const -> bool { return do_dups_ == 1; }
@@ -46,6 +47,7 @@ class run_mode {
   auto set_do_bisulfite(const int x) { if (x) do_bisulfite_ = x; }
   auto set_do_preseq(const int x) { if (x) do_preseq_ = x; }
   auto set_do_original_dups(const int x) { if (x) do_original_dups_ = x; }
+  auto set_do_quals(const int x) { if (x) do_quals_ = x; }
   //
   auto set_do_adap(const int x) { if (x) do_adap_ = x; }
   auto set_do_dups(const int x) { if (x) do_dups_ = x; }
@@ -74,6 +76,7 @@ class run_mode {
   static constexpr auto do_bisulfite_default = -1;      // OFF
   static constexpr auto do_preseq_default = -1;         // OFF
   static constexpr auto do_original_dups_default = -1;  // OFF
+  static constexpr auto do_quals_default = 1;           // ON
 
   static constexpr auto do_adap_default = 1;  // affects processing
   static constexpr auto do_dups_default = 1;  // affects processing
@@ -93,6 +96,7 @@ class run_mode {
   int do_bisulfite_{};
   int do_preseq_{};
   int do_original_dups_{};
+  int do_quals_{};
 
   int do_adap_{};
   int do_dups_{};
diff --git a/src/sam_file.cpp b/src/sam_file.cpp
index 4624c5c..f5c5d11 100644
--- a/src/sam_file.cpp
+++ b/src/sam_file.cpp
@@ -58,10 +58,11 @@ sam_file::shift_output_buffer() -> void {
 }
 
 auto
-sam_file::get_chunks(const std::int64_t n_chunks,  //
-                     const std::int32_t file_id,   //
-                     task_queue &tq,               //
+sam_file::get_chunks(std::int64_t n_chunks,
+                     const std::int32_t file_id,
+                     task_queue &tq,
                      std::atomic_int32_t &n_tasks) -> void {
+  static constexpr std::int64_t min_chunk_size = 64 * 1024;
   assert(n_chunks > 0);
   const auto beg_itr = std::begin(buffer);
   const auto end_itr = last;
@@ -78,7 +79,9 @@ sam_file::get_chunks(const std::int64_t n_chunks,  //
   };
   // clang-format on
   const std::int64_t n_bytes_available = std::distance(beg_itr, end_itr);
-  const auto [chunk_size, remainder] = std::div(n_bytes_available, n_chunks);
+  n_chunks = std::min((n_bytes_available + min_chunk_size - 1) / min_chunk_size,
+                      n_chunks);
+  auto [chunk_size, remainder] = std::div(n_bytes_available, n_chunks);
   assert(n_chunks > 0);
   auto start_itr = beg_itr;
   auto chunk_end = start_itr;
@@ -110,11 +113,12 @@ sam_file::load_next() -> void {
 }
 
 auto
-sam_file::make_tasks(const std::int64_t n_chunks,  //
-                     const std::int32_t file_id,   //
-                     task_queue &tq,               //
+sam_file::make_tasks(const std::int64_t n_threads,
+                     const std::int32_t file_id,
+                     task_queue &tq,
                      std::atomic_int32_t &n_tasks) -> void {
-  n_tasks = 1;  // for current task, which makes more tasks
+  static constexpr auto n_chunks_per_thread = 8;
+  const auto n_chunks = n_chunks_per_thread * n_threads;
   shift_output_buffer();
   load_next();
   get_chunks(n_chunks, file_id, tq, n_tasks);
diff --git a/src/sam_file.hpp b/src/sam_file.hpp
index a5c5358..fcc1621 100644
--- a/src/sam_file.hpp
+++ b/src/sam_file.hpp
@@ -15,7 +15,7 @@
 struct task_queue;
 
 class sam_file {
-  static constexpr auto min_buf_size = 64 * 1024;
+  static constexpr auto min_buf_size = 256 * 1024;
   std::vector<char> buffer;
   std::vector<char>::iterator cursor;
   std::vector<char>::iterator last;
@@ -33,12 +33,17 @@ class sam_file {
   ~sam_file() = default;
   // clang-format on
 
-  auto
-  make_tasks(const std::int64_t n_chunks,
+  friend auto
+  reset(sam_file &reads_file) -> void;
+
+  friend auto
+  make_tasks(sam_file &reads_file,
+             const std::int64_t n_threads,
              const std::int32_t file_id,
              task_queue &tq,
              std::atomic_int32_t &n_tasks) -> void;
 
+private:
   auto
   reset() -> void {
     buffer.clear();
@@ -47,13 +52,18 @@ class sam_file {
     last = cursor;
   }
 
-private:
   auto
   get_chunks(const std::int64_t n_chunks,
              const std::int32_t file_id,
              task_queue &tq,
              std::atomic_int32_t &n_tasks) -> void;
 
+  auto
+  make_tasks(const std::int64_t n_chunks,
+             const std::int32_t file_id,
+             task_queue &tq,
+             std::atomic_int32_t &n_tasks) -> void;
+
   auto
   shift_output_buffer() -> void;
 
@@ -70,9 +80,8 @@ make_tasks(sam_file &reads_file,
            const std::int32_t file_id,
            task_queue &tq,
            std::atomic_int32_t &n_tasks) -> void {
-  static constexpr auto n_chunks_per_thread = 8;
-  const auto n_chunks = n_chunks_per_thread * n_threads;
-  reads_file.make_tasks(n_chunks, file_id, tq, n_tasks);
+  n_tasks = 1;  // for current task, which makes more tasks
+  reads_file.make_tasks(n_threads, file_id, tq, n_tasks);
 }
 
 inline auto
diff --git a/src/sam_stdin.cpp b/src/sam_stdin.cpp
index 6da6f35..f4a881b 100644
--- a/src/sam_stdin.cpp
+++ b/src/sam_stdin.cpp
@@ -57,10 +57,11 @@ sam_stdin::skip_header() -> void {
 }
 
 auto
-sam_stdin::get_chunks(const std::int64_t n_chunks,
+sam_stdin::get_chunks(std::int64_t n_chunks,
                       const std::int32_t file_id,
                       task_queue &tq,
                       std::atomic_int32_t &n_tasks) -> void {
+  static constexpr std::int64_t min_chunk_size = 64 * 1024;
   assert(n_chunks > 0);
   const auto beg_itr = std::begin(buffer);
   const auto end_itr = last;
@@ -77,6 +78,8 @@ sam_stdin::get_chunks(const std::int64_t n_chunks,
   };
   // clang-format on
   const std::int64_t n_bytes_available = std::distance(beg_itr, end_itr);
+  n_chunks = std::min((n_bytes_available + min_chunk_size - 1) / min_chunk_size,
+                      n_chunks);
   const auto [chunk_size, remainder] = std::div(n_bytes_available, n_chunks);
   assert(n_chunks > 0);
   auto start_itr = beg_itr;
@@ -124,7 +127,7 @@ sam_stdin::load_next() -> void {
 validate_sam(const auto &buffer) -> bool {
   // SAM format has 11+ fields, tab separated and the docs give a regex for each
   // field. We are only checking the first char of the first field.
-  static constexpr auto n_bytes_to_validate = 16L * 1024;
+  static constexpr auto n_bytes_to_validate = 16 * 1024L;
   assert(std::size(buffer));
   const auto n_bytes = std::min(n_bytes_to_validate, std::ssize(buffer));
   const auto buf_end = std::cbegin(buffer) + n_bytes - (n_bytes >= 1L);
@@ -135,11 +138,12 @@ validate_sam(const auto &buffer) -> bool {
 }
 
 auto
-sam_stdin::make_tasks(const std::int64_t n_chunks,
+sam_stdin::make_tasks(const std::int64_t n_threads,
                       const std::int32_t file_id,
                       task_queue &tq,
                       std::atomic_int32_t &n_tasks) -> void {
-  n_tasks = 1;  // for current task, which makes more tasks
+  static constexpr auto n_chunks_per_thread = 8;
+  const auto n_chunks = n_chunks_per_thread * n_threads;
   shift_output_buffer();
   load_next();
   if (!validate_sam(buffer))
diff --git a/src/sam_stdin.hpp b/src/sam_stdin.hpp
index c9f0bdd..1ee7275 100644
--- a/src/sam_stdin.hpp
+++ b/src/sam_stdin.hpp
@@ -14,7 +14,7 @@
 struct task_queue;
 
 class sam_stdin {
-  static constexpr auto min_buf_size = 64 * 1024;
+  static constexpr auto min_buf_size = 256 * 1024;
   std::vector<char> buffer;
   std::vector<char>::iterator cursor;
   std::vector<char>::iterator last;
@@ -32,8 +32,25 @@ class sam_stdin {
   ~sam_stdin() = default;
   // clang-format on
 
+  friend auto
+  reset(sam_stdin &reads_file) -> void;
+
+  friend auto
+  make_tasks(sam_stdin &reads_file,
+             const std::int64_t n_threads,
+             const std::int32_t file_id,
+             task_queue &tq,
+             std::atomic_int32_t &n_tasks) -> void;
+
+private:
   auto
-  make_tasks(const std::int64_t n_chunks,
+  get_chunks(std::int64_t n_chunks,
+             const std::int32_t file_id,
+             task_queue &tq,
+             std::atomic_int32_t &n_tasks) -> void;
+
+  auto
+  make_tasks(const std::int64_t n_threads,
              const std::int32_t file_id,
              task_queue &tq,
              std::atomic_int32_t &n_tasks) -> void;
@@ -46,13 +63,6 @@ class sam_stdin {
     last = cursor;
   }
 
-private:
-  auto
-  get_chunks(const std::int64_t n_chunks,
-             const std::int32_t file_id,
-             task_queue &tq,
-             std::atomic_int32_t &n_tasks) -> void;
-
   auto
   shift_output_buffer() -> void;
 
@@ -73,9 +83,8 @@ make_tasks(sam_stdin &reads_file,
            const std::int32_t file_id,
            task_queue &tq,
            std::atomic_int32_t &n_tasks) -> void {
-  static constexpr auto n_chunks_per_thread = 8;
-  const auto n_chunks = n_chunks_per_thread * n_threads;
-  reads_file.make_tasks(n_chunks, file_id, tq, n_tasks);
+  n_tasks = 1;  // for current task, which makes more tasks
+  reads_file.make_tasks(n_threads, file_id, tq, n_tasks);
 }
 
 inline auto
diff --git a/src/samrec.cpp b/src/samrec.cpp
index 54bab1e..98b4bb4 100644
--- a/src/samrec.cpp
+++ b/src/samrec.cpp
@@ -92,7 +92,7 @@ samrec::get_next(samrec::pos_t &cursor,
   if (itr++ == end_itr)
     return false;
 
-  if (qual_len != seq_len)
+  if (qual_len != seq_len && qual_len != 1)
     throw std::runtime_error("quality scores invalid for record: " +
                              std::string(name_itr, name_itr + name_len));
 
@@ -108,13 +108,20 @@ samrec::get_next(samrec::pos_t &cursor,
   if (bam_is_rev(flag)) {
     std::reverse_copy(seq_itr, seq_itr + seq_len, out_seq_itr);
     std::transform(out_seq_itr, out_seq_itr + seq_len, out_seq_itr, complem);
-    std::reverse_copy(qual_itr, qual_itr + seq_len, out_qual_itr);
+    if (*qual_itr != qual_missing_code)
+      std::reverse_copy(qual_itr, qual_itr + seq_len, out_qual_itr);
+    else
+      *out_qual_itr = qual_missing_code;
   }
   else {
     std::copy(seq_itr, seq_itr + seq_len, out_seq_itr);
-    std::copy(qual_itr, qual_itr + seq_len, out_qual_itr);
+    if (*qual_itr != qual_missing_code)
+      std::copy(qual_itr, qual_itr + seq_len, out_qual_itr);
+    else
+      *out_qual_itr = qual_missing_code;
   }
-  assert(std::all_of(out_qual_itr, out_qual_itr + seq_len, [](const auto q) {
+
+  assert(std::all_of(get_qual(rec), get_qual_end(rec), [](const auto q) {
     return q >= 0 && q <= falco::max_qual_val;
   }));
 
diff --git a/src/samrec.hpp b/src/samrec.hpp
index 3e96f8d..ce48f0b 100644
--- a/src/samrec.hpp
+++ b/src/samrec.hpp
@@ -27,6 +27,7 @@ class samrec {
   friend constexpr auto get_seq(const samrec &);
   friend constexpr auto get_seq_end(const samrec &);
   friend constexpr auto get_seq_size(const samrec &);
+  friend constexpr auto has_qual(const samrec &);
   friend constexpr auto get_qual(const samrec &);
   friend constexpr auto get_qual_end(const samrec &);
   friend constexpr auto get_qual_size(const samrec &);
@@ -49,6 +50,11 @@ class samrec {
 
   [[nodiscard]] static auto
   find_end_pos(pos_t itr, const pos_t end) -> pos_t;
+
+  [[nodiscard]] static consteval auto
+  get_qual_missing_code() -> char {
+    return qual_missing_code;
+  }
 };
 
 [[nodiscard]] inline constexpr auto
@@ -83,13 +89,18 @@ get_qual(const samrec &rec) {
 }
 
 [[nodiscard]] inline constexpr auto
-get_qual_end(const samrec &rec) {
-  return get_qual(rec) + rec.seq_len;  // NOLINT
+has_qual(const samrec &rec) {
+  return *get_qual(rec) != samrec::qual_missing_code;
 }
 
 [[nodiscard]] inline constexpr auto
 get_qual_size(const samrec &rec) {
-  return get_seq_size(rec);
+  return has_qual(rec) ? rec.seq_len : 0U;
+}
+
+[[nodiscard]] inline constexpr auto
+get_qual_end(const samrec &rec) {
+  return get_qual(rec) + get_qual_size(rec);
 }
 
 struct sam_task_t {
diff --git a/src/tile_processor.cpp b/src/tile_processor.cpp
index 5fe1d93..e2386ac 100644
--- a/src/tile_processor.cpp
+++ b/src/tile_processor.cpp
@@ -205,6 +205,7 @@ tile_processor::add_and_consume(
     else
       quals.emplace(rhs_tile_id, std::move(rhs_qual));
   }
+  max_read_len = std::max(max_read_len, rhs.max_read_len);
   rhs.release();
 }