DAY-EC activated.
DAY-EC activated.
repo=/fsx/analysis_results/ubuntu/hybonly_ilmn_kitchensink_mounted_20260606T053415Z/daylily-omics-analysis
head=d02a1bd825c8c13bac3751d26854b9834c5c3a53
branch=jem-dev
status:
## jem-dev...origin/jem-dev
 M workflow/rules/contam_identity.smk
 M workflow/rules/go_left.smk
 M workflow/rules/mosdepth.smk
 M workflow/rules/rtg_vcfeval.smk
 M workflow/rules/vep.smk
?? pipeline_details.md
?? pipeline_workflow_checkpoint_20260606T055941Z.mmd
?? pipeline_workflow_checkpoint_20260606T055941Z.pdf
?? pipeline_workflow_checkpoint_20260606T061441Z.mmd
?? pipeline_workflow_checkpoint_20260606T061441Z.pdf
?? pipeline_workflow_checkpoint_20260606T062942Z.mmd
?? pipeline_workflow_checkpoint_20260606T062942Z.pdf
?? pipeline_workflow_checkpoint_20260606T064443Z.mmd
?? pipeline_workflow_checkpoint_20260606T064443Z.pdf
?? pipeline_workflow_checkpoint_20260606T065944Z.mmd
?? pipeline_workflow_checkpoint_20260606T065944Z.pdf
?? pipeline_workflow_checkpoint_20260606T071445Z.mmd
?? pipeline_workflow_checkpoint_20260606T071445Z.pdf
?? pipeline_workflow_checkpoint_20260606T073009Z.mmd
?? pipeline_workflow_checkpoint_20260606T073009Z.pdf
?? pipeline_workflow_checkpoint_20260606T074509Z.mmd
?? pipeline_workflow_checkpoint_20260606T074509Z.pdf
?? pipeline_workflow_checkpoint_20260606T080010Z.mmd
?? pipeline_workflow_checkpoint_20260606T080010Z.pdf
?? pipeline_workflow_checkpoint_20260606T081511Z.mmd
?? pipeline_workflow_checkpoint_20260606T081511Z.pdf
?? pipeline_workflow_checkpoint_20260606T083012Z.mmd
?? pipeline_workflow_checkpoint_20260606T083012Z.pdf
?? pipeline_workflow_checkpoint_20260606T084513Z.mmd
?? pipeline_workflow_checkpoint_20260606T084513Z.pdf
?? pipeline_workflow_checkpoint_20260606T090014Z.mmd
?? pipeline_workflow_checkpoint_20260606T090014Z.pdf
?? pipeline_workflow_planned.mmd
?? pipeline_workflow_planned.pdf
?? sbatch_errs.log
diff_stat:
 workflow/rules/contam_identity.smk | 95 ++++++++++++++++++++++++++++++++++----
 workflow/rules/go_left.smk         | 12 ++++-
 workflow/rules/mosdepth.smk        | 15 ++++--
 workflow/rules/rtg_vcfeval.smk     |  5 +-
 workflow/rules/vep.smk             |  6 ++-
 5 files changed, 116 insertions(+), 17 deletions(-)
diff_target:
diff --git a/workflow/rules/contam_identity.smk b/workflow/rules/contam_identity.smk
index fbb3d26..54e4d9f 100644
--- a/workflow/rules/contam_identity.smk
+++ b/workflow/rules/contam_identity.smk
@@ -328,13 +328,51 @@ rule haplocheck_vcf_contam_identity:
             rm -rf "$result_dir"
         }}
         trap cleanup EXIT
-        {params.command:q} --out "$result_prefix" --raw {input.vcf:q} > {log:q} 2>&1
-        test -s "$result_dir/contamination.txt"
-        test -s "$result_dir/contamination.raw.txt"
-        test -s "$result_dir/contamination.html"
-        cp "$result_dir/contamination.txt" {output.contamination:q}
-        cp "$result_dir/contamination.raw.txt" {output.raw:q}
-        cp "$result_dir/contamination.html" {output.html:q}
+        set +o pipefail
+        if gzip -cd {input.vcf:q} | grep -m 1 -q -v '^#'; then
+            has_variants=true
+        else
+            has_variants=false
+        fi
+        set -o pipefail
+        if [[ "$has_variants" == "true" ]]; then
+            set +e
+            {params.command:q} --out "$result_prefix" --raw {input.vcf:q} > {log:q} 2>&1
+            haplocheck_rc=$?
+            set -e
+            if grep -q 'outside the range.*rCRS only' {log:q}; then
+                printf 'SampleID	Contamination Status	Contamination Level	Distance	Sample Coverage	Major Haplogroup	Minor Haplogroup
+%s	UNSUPPORTED_REFERENCE	0		0		
+' {wildcards.sample:q} > {output.contamination:q}
+                printf 'SampleID	Contamination Status	Contamination Level	Distance	Sample Coverage	Major Haplogroup	Minor Haplogroup
+%s	UNSUPPORTED_REFERENCE	0		0		
+' {wildcards.sample:q} > {output.raw:q}
+                printf '<html><body>UNSUPPORTED_REFERENCE</body></html>
+' > {output.html:q}
+                printf 'UNSUPPORTED_REFERENCE: haplocheck skipped because the input VCF is not restricted to rCRS positions.
+' >> {log:q}
+            elif [[ "$haplocheck_rc" -eq 0 ]]; then
+                test -s "$result_dir/contamination.txt"
+                test -s "$result_dir/contamination.raw.txt"
+                test -s "$result_dir/contamination.html"
+                cp "$result_dir/contamination.txt" {output.contamination:q}
+                cp "$result_dir/contamination.raw.txt" {output.raw:q}
+                cp "$result_dir/contamination.html" {output.html:q}
+            else
+                exit "$haplocheck_rc"
+            fi
+        else
+            printf 'SampleID	Contamination Status	Contamination Level	Distance	Sample Coverage	Major Haplogroup	Minor Haplogroup
+%s	NO_VARIANTS	0		0		
+' {wildcards.sample:q} > {output.contamination:q}
+            printf 'SampleID	Contamination Status	Contamination Level	Distance	Sample Coverage	Major Haplogroup	Minor Haplogroup
+%s	NO_VARIANTS	0		0		
+' {wildcards.sample:q} > {output.raw:q}
+            printf '<html><body>NO_VARIANTS</body></html>
+' > {output.html:q}
+            printf 'NO_VARIANTS: haplocheck skipped because the input VCF has no variant records.
+' > {log:q}
+        fi
         """
 
 
@@ -371,9 +409,46 @@ rule read_haps_contam_identity:
         set -euo pipefail
         test {params.sample_ok:q} = ok
         mkdir -p "$(dirname {output.txt:q})" "$(dirname {log:q})"
-        command -v {params.command:q} > /dev/null
-        test -s {params.reliable_snp_file:q}
-        {params.command:q} {params.extra_args} -fa {params.ref:q} {input.bam:q} {params.reliable_snp_file:q} {input.vcf:q} > {output.txt:q} 2> {log:q}
+        set +o pipefail
+        if gzip -cd {input.vcf:q} | grep -m 1 -q -v '^#'; then
+            has_variants=true
+        else
+            has_variants=false
+        fi
+        set -o pipefail
+        if [[ "$has_variants" == "true" ]]; then
+            if ! command -v {params.command:q} > /dev/null; then
+                printf 'SNP_PAIRS ERROR_PAIRS DOUBLE_ERROR_PAIR_COUNT DOUBLE_ERROR_FRACTION REL_ERROR_FRACTION NONSENSE_FRACTION PASS_FAIL REASON
+0 0 0 0 0 0 NO_DATA READ_HAPS_UNAVAILABLE
+' > {output.txt:q}
+                printf 'READ_HAPS_UNAVAILABLE: read_haps command is unavailable: %s
+' {params.command:q} > {log:q}
+            elif [[ ! -s {params.reliable_snp_file:q} ]]; then
+                printf 'SNP_PAIRS ERROR_PAIRS DOUBLE_ERROR_PAIR_COUNT DOUBLE_ERROR_FRACTION REL_ERROR_FRACTION NONSENSE_FRACTION PASS_FAIL REASON
+0 0 0 0 0 0 NO_DATA READ_HAPS_MARKERS_UNAVAILABLE
+' > {output.txt:q}
+                printf 'READ_HAPS_MARKERS_UNAVAILABLE: read_haps marker file is missing or empty: %s
+' {params.reliable_snp_file:q} > {log:q}
+            else
+                set +e
+                {params.command:q} {params.extra_args} -fa {params.ref:q} {input.bam:q} {params.reliable_snp_file:q} {input.vcf:q} > {output.txt:q} 2> {log:q}
+                read_haps_rc=$?
+                set -e
+                if [[ "$read_haps_rc" != "0" ]] || [[ ! -s {output.txt:q} ]] || ! grep -q 'PASS_FAIL' {output.txt:q} || ! grep -q 'REASON' {output.txt:q}; then
+                printf 'SNP_PAIRS ERROR_PAIRS DOUBLE_ERROR_PAIR_COUNT DOUBLE_ERROR_FRACTION REL_ERROR_FRACTION NONSENSE_FRACTION PASS_FAIL REASON
+0 0 0 0 0 0 NO_DATA READ_HAPS_FAILED
+' > {output.txt:q}
+                    printf 'READ_HAPS_FAILED: read_haps exited with status %s or wrote no usable QC table.
+' "$read_haps_rc" >> {log:q}
+                fi
+            fi
+        else
+            printf 'SNP_PAIRS ERROR_PAIRS DOUBLE_ERROR_PAIR_COUNT DOUBLE_ERROR_FRACTION REL_ERROR_FRACTION NONSENSE_FRACTION PASS_FAIL REASON
+0 0 0 0 0 0 NO_DATA NO_VARIANTS
+' > {output.txt:q}
+            printf 'NO_VARIANTS: read_haps skipped because the input VCF has no variant records.
+' > {log:q}
+        fi
         test -s {output.txt:q}
         grep -q 'PASS_FAIL' {output.txt:q}
         grep -q 'REASON' {output.txt:q}
