Add read counting to demux_fastqs task and update workflows

- Add samtools view -c to count reads in demux_fastqs task
- Add read_counts Array[Int] output to demux_fastqs task
- Add read_counts_raw output to load_illumina_fastqs workflow
- Update load_illumina_fastqs_deplete to use demux read counts
  instead of spikein counts (which default to 0 when no spikein_db)
- Move test input/output JSONs to test/input/WDL/ with symlinks
  to miniwdl-local/ and cromwell-local/
- Add read_counts_raw validation to test output JSONs
- Remove outdated comments about viral-core bug #127 (now fixed)

🤖 Generated with [Claude Code](https://claude.com/claude-code)

Co-Authored-By: Claude <noreply@anthropic.com>

Author: dpark01

pipes/WDL/tasks/tasks_demux.wdl

@@ -893,9 +893,6 @@ task get_illumina_run_metadata {
 
     illumina.py --version | tee VERSION
 
-    # NOTE: This command currently fails in viral-core:2.5.1 due to bug #127
-    # See: https://github.com/broadinstitute/viral-core/issues/127
-    # Will work once viral-core is updated
     illumina.py illumina_metadata \
       --samplesheet ~{samplesheet} \
       --runinfo ~{runinfo_xml} \
@@ -976,10 +973,6 @@ task demux_fastqs {
 
     illumina.py --version | tee VERSION
 
-    # NOTE: This command currently fails in viral-core:2.5.1 due to bug #127
-    # See: https://github.com/broadinstitute/viral-core/issues/127
-    # Also note: run_date and flowcell_id are extracted from RunInfo.xml internally
-    # and cannot be overridden (feature request in same issue)
     illumina.py splitcode_demux_fastqs \
       --fastq_r1 ~{fastq_r1} \
       ~{'--fastq_r2 ' + fastq_r2} \
@@ -992,9 +985,16 @@ task demux_fastqs {
     # Count output BAMs
     ls -lh *.bam || echo "No BAM files found"
 
-    # Run FastQC on output BAMs (optional, if BAMs were created)
+    # Initialize read counts file (empty in case no BAMs are created)
+    touch read_counts.txt
+
+    # Run FastQC and count reads for output BAMs
     if ls *.bam 1> /dev/null 2>&1; then
       for bam in *.bam; do
+        # Count reads in BAM and append to read counts file
+        samtools view -c "$bam" >> read_counts.txt
+
+        # Run FastQC
         reports.py fastqc \
           "$bam" \
           "${bam%.bam}_fastqc.html" \
@@ -1007,6 +1007,7 @@ task demux_fastqs {
 
   output {
     Array[File] output_bams      = glob("*.bam")
+    Array[Int]  read_counts      = read_lines("read_counts.txt")
     Array[File] fastqc_html      = glob("*_fastqc.html")
     Array[File] fastqc_zip       = glob("*_fastqc.zip")
     String      viralngs_version = read_string("VERSION")

pipes/WDL/workflows/load_illumina_fastqs.wdl

@@ -55,6 +55,7 @@ workflow load_illumina_fastqs {
   output {
     # BAM outputs (flattened)
     Array[File] raw_reads_unaligned_bams = flatten(demux_fastqs.output_bams)
+    Array[Int]  read_counts_raw          = flatten(demux_fastqs.read_counts)
 
     # QC outputs (flattened)
     Array[File] fastqc_html = flatten(demux_fastqs.fastqc_html)

pipes/WDL/workflows/load_illumina_fastqs_deplete.wdl

@@ -86,9 +86,13 @@ workflow load_illumina_fastqs_deplete {
   }
 
   Array[File] raw_bams = flatten(demux_fastqs.output_bams)
+  Array[Int]  raw_read_counts = flatten(demux_fastqs.read_counts)
 
   # Step 6: Spike-in counting and depletion for all BAMs
-  scatter (raw_reads in raw_bams) {
+  scatter (idx in range(length(raw_bams))) {
+    File raw_reads = raw_bams[idx]
+    Int  raw_count = raw_read_counts[idx]
+
     # Spike-in counting (if spikein_db provided)
     if (defined(spikein_db)) {
       call reports.align_and_count as spikein {
@@ -97,7 +101,6 @@ workflow load_illumina_fastqs_deplete {
           ref_db = select_first([spikein_db])
       }
     }
-    Int reads_total_count = select_first([spikein.reads_total, 0])
 
     # Depletion (if any depletion db provided)
     if (length(flatten(select_all([bmtaggerDbs, blastDbs, bwaDbs]))) > 0) {
@@ -110,7 +113,7 @@ workflow load_illumina_fastqs_deplete {
       }
     }
 
-    Int  read_count_post_depletion = select_first([deplete.depletion_read_count_post, reads_total_count])
+    Int  read_count_post_depletion = select_first([deplete.depletion_read_count_post, raw_count])
     File cleaned_bam = select_first([deplete.cleaned_bam, raw_reads])
 
     if (read_count_post_depletion >= min_reads_per_bam) {
@@ -120,7 +123,7 @@ workflow load_illumina_fastqs_deplete {
       File empty_bam = raw_reads
     }
 
-    Pair[String,Int] count_raw = (basename(raw_reads, '.bam'), reads_total_count)
+    Pair[String,Int] count_raw = (basename(raw_reads, '.bam'), raw_count)
     Pair[String,Int] count_cleaned = (basename(raw_reads, '.bam'), read_count_post_depletion)
   }
 
@@ -212,7 +215,7 @@ workflow load_illumina_fastqs_deplete {
   output {
     # Raw BAM outputs
     Array[File] raw_reads_unaligned_bams = raw_bams
-    Array[Int]  read_counts_raw          = reads_total_count
+    Array[Int]  read_counts_raw          = raw_read_counts
 
     # Metadata outputs (with defaults applied)
     Map[String,Map[String,String]] meta_by_sample   = read_json(meta_default_sample.out_json)

test/input/WDL/cromwell-local/test_inputs-load_illumina_fastqs-local.json

@@ -0,0 +1 @@
+../test_inputs-load_illumina_fastqs-local.json

test/input/WDL/cromwell-local/test_inputs-load_illumina_fastqs_deplete-local.json

@@ -0,0 +1 @@
+../test_inputs-load_illumina_fastqs_deplete-local.json

test/input/WDL/cromwell-local/test_outputs-load_illumina_fastqs-local.json

@@ -0,0 +1 @@
+../test_outputs-load_illumina_fastqs-local.json

test/input/WDL/cromwell-local/test_outputs-load_illumina_fastqs_deplete-local.json

@@ -0,0 +1 @@
+../test_outputs-load_illumina_fastqs_deplete-local.json

test/input/WDL/miniwdl-local/test_inputs-load_illumina_fastqs-local.json

@@ -0,0 +1 @@
+../test_inputs-load_illumina_fastqs-local.json