# shellcheck shell=bash
# common.sh: sourced by every script in scripts/, once SAMPLE and GENOME_DIR
# are read from the arguments and the environment:
#   # shellcheck source=lib/common.sh
#   . "$(dirname "$0")/lib/common.sh"
#   validate_sample "$SAMPLE"
#
# It gives each script the same:
#   versions.env   every image tag and coupled data version, sourced with no
#                  fallback: a missing file or a typo stops the script.
#   umask 077      files the host side writes are readable by their owner only.
#   GENOME_DIR     as exported by the caller (not required here: setup.sh and
#                  validate-setup.sh run without it).
#   REF_FASTA      host path of the GRCh38 FASTA. Default: the NCBI GRCh38
#                  no-ALT analysis set setup.sh installs. Override from the
#                  environment, absolute or relative to GENOME_DIR; it must lie
#                  inside GENOME_DIR, the only directory containers see.
#                  REF_FASTA_C is the same file inside a container, REF_DICT
#                  its GATK sequence dictionary.
#   THREADS        CPU budget of one step. Default 8; a script that wants another
#                  default sets THREADS=${THREADS:-N} before sourcing this file.
#   CONTAINER_ENGINE  docker by default.
# and the helpers validate_sample, require_image, cpath, run_in, have_output,
# wrote_vcf, atomic_out, fetch, data_file / install_data_file, install_vep_cache,
# lock_acquire / lock_release, pipeline_images, pinned_tag / tag_pinned_image /
# pinned_tag_present, image_is_dangling,
# scatter_beds and run_parallel, FASTQ_RECORDS_AWK / fastp_total_reads /
# read_count_mode / read_count_check,
# described where they are
# defined. Keep it bash 3.2 compatible: macOS runs setup.sh with /bin/bash.

PGP_ROOT=$(cd "$(dirname "${BASH_SOURCE[0]}")/../.." && pwd)
# shellcheck source=../versions.env
. "${PGP_ROOT}/versions.env"

umask 077
CONTAINER_ENGINE=${CONTAINER_ENGINE:-docker}
THREADS=${THREADS:-8}

# validate_sample NAME: exit 2 unless NAME is a plain sample name. Scripts put
# it into container paths and `bash -c` bodies, so only letters, digits, '.',
# '_' and '-' are allowed (the same rule as main.nf), and '.' and '..' are not
# names at all. Nor is a name that starts with '-': CPSR, bin/pgx_parse.py and
# bin/collect_summary.py take it as an argument value, and argparse reads '-x'
# as an option. CPSR's 3 to 40 characters are bin/cpsr_sample_id's job.
validate_sample() {
  local s=${1:-}
  if [[ ! "$s" =~ ^[A-Za-z0-9._-]+$ ]] || [ "$s" = "." ] || [ "$s" = ".." ]; then
    echo "ERROR: invalid sample name '${s}'. Use only letters, digits, '.', '_' and '-' (and not '.' or '..')." >&2
    exit 2
  fi
  if [[ "$s" == -* ]]; then
    echo "ERROR: invalid sample name '${s}': it starts with '-', which tools read as an option. Start it with a letter or digit." >&2
    exit 2
  fi
}

# require_image VAR...: exit 1 unless each VAR is set by versions.env.
require_image() {
  local v
  for v in "$@"; do
    if [ -z "${!v:-}" ]; then
      echo "ERROR: ${v} is not set; add it to versions.env." >&2
      exit 1
    fi
  done
}

# cpath HOST_PATH: the path a container sees for a host path inside GENOME_DIR.
cpath() {
  case "$1" in
    "${GENOME_DIR:?}") printf '/genome' ;;
    "${GENOME_DIR}"/*) printf '/genome%s' "${1#"$GENOME_DIR"}" ;;
    *) echo "ERROR: ${1} is outside GENOME_DIR (${GENOME_DIR}); containers only see GENOME_DIR." >&2
       return 1 ;;
  esac
}

# shellcheck disable=SC2034  # REF_DICT and REF_FASTA_C are read by the scripts
if [ -n "${GENOME_DIR:-}" ]; then
  GENOME_DIR=${GENOME_DIR%/}
  # A new file name on purpose: a reference with ALT contigs left from an
  # older version is never picked up by accident (docs/realignment.md).
  REF_FASTA=${REF_FASTA:-${GENOME_DIR}/reference/GRCh38_no_alt_analysis_set.fasta}
  case "$REF_FASTA" in
    /*) ;;
    *) REF_FASTA="${GENOME_DIR}/${REF_FASTA}" ;;
  esac
  # GATK names the dictionary after the FASTA without its extensions:
  # ref.fasta and ref.fa.gz both give ref.dict.
  REF_DICT="${REF_FASTA%.gz}"
  REF_DICT="${REF_DICT%.*}.dict"
  REF_FASTA_C=$(cpath "$REF_FASTA") || exit 2
fi

# run_in [--net] [--root] [--rw DIR]... [docker run options] IMAGE [COMMAND...]
#
# `docker run --rm` with these defaults:
#   --network none            the analysis steps need no network; a tool that
#                             reaches out at run time fails loudly instead.
#   GENOME_DIR at /genome:ro  references, ClinVar, caches and bundles are
#                             read-only, and only
#   GENOME_DIR/SAMPLE         (when SAMPLE is set) is writable, at /genome/SAMPLE.
#   --user UID:GID, HOME=/tmp outputs belong to the caller, not to root.
# Opt-outs, each named at the call site with the reason:
#   --net      the step downloads something (a database, a pip package).
#   --root     the image cannot run as an unprivileged user: the container runs
#              as root (--user 0:0), also when the image names another user.
#   --rw DIR   DIR (inside GENOME_DIR) is writable too, e.g. a shared index.
# Everything after the opt-outs goes to `docker run` unchanged.
# When an unprivileged run fails and the sample directory holds a path the
# caller does not own (versions before this library ran every container as
# root), run_in says how to take the directory back.
run_in() {
  local isolate=true root=false d c rc=0 other
  local -a extra=()
  while [ $# -gt 0 ]; do
    case "$1" in
      --net) isolate=false; shift ;;
      --root) root=true; shift ;;
      --rw)
        d=${2:?run_in --rw needs a directory}
        c=$(cpath "$d") || exit 2
        mkdir -p "$d"
        extra+=(-v "${d}:${c}")
        shift 2 ;;
      *) break ;;
    esac
  done
  local -a args=(run --rm)
  if $isolate; then args+=(--network none); fi
  # --root is explicit: the container runs as root (--user 0:0) even when
  # the image names another user. It is for an image that has not been shown
  # to run unprivileged, for example PCGR in scripts/17-cpsr.sh.
  if $root; then
    args+=(--user 0:0)
  else
    args+=(--user "$(id -u):$(id -g)" -e HOME=/tmp)
  fi
  args+=(-v "${GENOME_DIR:?}:/genome:ro")
  if [ -n "${SAMPLE:-}" ]; then
    mkdir -p "${GENOME_DIR}/${SAMPLE}"
    args+=(-v "${GENOME_DIR}/${SAMPLE}:/genome/${SAMPLE}")
  fi
  "$CONTAINER_ENGINE" "${args[@]}" ${extra[@]+"${extra[@]}"} "$@" || rc=$?
  if [ "$rc" -ne 0 ] && ! $root && [ -n "${SAMPLE:-}" ]; then
    # A root-owned *.lock (left by a Nextflow CRAM_ARCHIVE task under
    # --user root) is opened read-only and blocks nothing: not worth the hint.
    other=$(find "${GENOME_DIR}/${SAMPLE}" ! -user "$(id -u)" ! -name '*.lock' -print 2>/dev/null | head -n 1)
    if [ -n "$other" ]; then
      echo "NOTE: ${other} is not yours. Older versions of this pipeline ran every container as root;" >&2
      echo "  steps now run as you and cannot write over such files. If the error above is" >&2
      echo "  'Permission denied', take the sample directory back and run the step again:" >&2
      echo "    sudo chown -R \"$(id -u):$(id -g)\" \"${GENOME_DIR}/${SAMPLE}\"" >&2
    fi
  fi
  return "$rc"
}

# --- The alignment read count (steps 02 and 02a, validate-setup.sh) ----------
# minimap2 and bwa-mem2 write every read they are given, mapped or not, as one
# primary record (neither 0x100 secondary nor 0x800 supplementary), and samtools
# fixmate, sort and markdup drop none. So a finished BAM holds as many primary
# records as R1 and R2 together hold FASTQ records: reads, not pairs, which is
# also the unit of fastp's summary.after_filtering.total_reads. Both aligners
# stop at the end of the shorter file of a pair with only a warning, so a cut
# R1 or R2 loses reads without an error. The ALIGN_MINIMAP2 module applies the
# same rule. READ_COUNT_CHECK=warn turns a difference into a warning.

# FASTQ_RECORDS_AWK: an awk program that prints the number of FASTQ records in
# its decompressed input, or exits 3 when the line count is not a multiple of
# 4 (a record cut short). %.0f, not %d: a 30x file has more than 2^31 lines.
# shellcheck disable=SC2034  # read by the step scripts
FASTQ_RECORDS_AWK='END { if (NR % 4) { printf "%.0f lines, not a multiple of 4: a record is cut short\n", NR > "/dev/stderr"; exit 3 } printf "%.0f\n", NR / 4 }'

# fastp_total_reads JSON: summary.after_filtering.total_reads of a fastp
# report (the reads fastp wrote, R1 plus R2 records). Exit status 1 when the
# report has none. The summary block comes first in fastp's report, so the
# first "after_filtering" object is the summary's.
fastp_total_reads() {
  tr -d ' \t\r\n' < "$1" | awk '
    { k = "\"after_filtering\":{\"total_reads\":"
      i = index($0, k)
      if (i) { v = substr($0, i + length(k)); sub(/[^0-9].*/, "", v); if (v != "") { print v; found = 1 } }
      exit }
    END { exit !found }'
}

# read_count_mode: prints READ_COUNT_CHECK (fail, the default, or warn), or
# says what is wrong and returns 2. A step calls it before it starts work.
read_count_mode() {
  case "${READ_COUNT_CHECK:-fail}" in
    fail|warn) echo "${READ_COUNT_CHECK:-fail}" ;;
    *) echo "ERROR: READ_COUNT_CHECK must be 'fail' or 'warn', got '${READ_COUNT_CHECK}'" >&2; return 2 ;;
  esac
}

# read_count_check PRIMARY SOURCE TABLE COUNT...: compare the primary records
# of a BAM with the reads that went into the aligner, the sum of the COUNTs
# (one per FASTQ, or fastp's total). SOURCE names where the COUNTs come from.
# Returns 0 when they are equal; 1 when they differ (what is missing goes to
# stderr); 3 when they differ and READ_COUNT_CHECK=warn (the same text, as a
# warning); 2 when a count is not a whole number or READ_COUNT_CHECK is bad.
# TABLE, unless empty, gets the key/value table the reports read: input_reads,
# input_source, primary_records, missing_reads (input minus primary) and
# read_count_check (ok or mismatch).
read_count_check() {
  local primary=$1 src=$2 table=$3 mode c total=0 missing verdict=ok level pct
  shift 3
  mode=$(read_count_mode) || return 2
  for c in "$primary" "$@"; do
    if [[ ! "$c" =~ ^(0|[1-9][0-9]*)$ ]]; then
      echo "ERROR: could not count the reads: '${primary}' primary records in the BAM, '$*' reads in (${src})" >&2
      return 2
    fi
  done
  [ "$#" -gt 0 ] || { echo "ERROR: read_count_check got no input count" >&2; return 2; }
  for c in "$@"; do total=$((total + c)); done
  missing=$((total - primary))
  [ "$missing" -eq 0 ] || verdict=mismatch
  if [ -n "$table" ]; then
    printf 'key\tvalue\ninput_reads\t%s\ninput_source\t%s\nprimary_records\t%s\nmissing_reads\t%s\nread_count_check\t%s\n' \
      "$total" "$src" "$primary" "$missing" "$verdict" > "$table"
  fi
  if [ "$verdict" = ok ]; then
    echo "Read count: the BAM holds ${primary} primary records, the ${total} reads that went in (${src})."
    return 0
  fi
  level=ERROR
  [ "$mode" = fail ] || level=WARNING
  if [ "$missing" -gt 0 ]; then
    pct=$(awk -v m="$missing" -v n="$total" 'BEGIN { printf "%.1f", 100 * m / n }')
    echo "${level}: the BAM holds ${primary} primary records, but ${total} reads went in (${src}): ${missing} reads (${pct}%) are missing." >&2
  else
    echo "${level}: the BAM holds ${primary} primary records, more than the ${total} reads that went in (${src}): it was not made from these reads." >&2
  fi
  echo "  The aligner writes every read it is given, mapped or not, once, so the two counts must be equal." >&2
  echo "  Depth, copy number (CYP2D6), CNV, telomere content and every call that depends on coverage would be off." >&2
  echo "  An R1 and R2 with different numbers of records is one cause: minimap2 and bwa-mem2 stop at the end" >&2
  echo "  of the shorter file with only a warning. See docs/troubleshooting.md, \"The BAM lost reads\"." >&2
  if [ "$mode" = warn ]; then
    echo "  READ_COUNT_CHECK=warn: going on." >&2
    return 3
  fi
  echo "  For an input you know is odd, READ_COUNT_CHECK=warn goes on anyway." >&2
  return 1
}

# have_output FILE...: true when every FILE is a finished output, the test a
# step uses before it skips work. Every FILE must be non-empty, and by name:
#   .vcf                  starts with ##fileformat=VCF and has a #CHROM line;
#   .vcf.gz, .bcf, .bam   ends with the BGZF end-of-file block, which htslib
#                         and htsjdk write only when they close the file, and
#                         starts with the format's own header.
# So an empty file left by a failed redirect, or one cut short by a killed
# run, is redone instead of trusted.
have_output() {
  local f first
  for f in "$@"; do
    [ -s "$f" ] || return 1
    case "$f" in
      *.vcf)
        first=$(head -n 1 "$f")
        [[ "$first" == "##fileformat=VCF"* ]] || return 1
        grep -q '^#CHROM' "$f" || return 1 ;;
      *.vcf.gz|*.bcf|*.bam)
        _bgzf_complete "$f" || return 1
        first=$(gzip -cd "$f" 2>/dev/null | head -c 16 | tr -d '\000') || true
        case "$f" in
          *.vcf.gz) [[ "$first" == "##fileformat=VCF"* ]] || return 1 ;;
          *.bcf) [[ "$first" == BCF* ]] || return 1 ;;
          *.bam) [[ "$first" == BAM* ]] || return 1 ;;
        esac ;;
    esac
  done
}

# wrote_vcf FILE: true when FILE is non-empty and starts with a VCF header,
# plain or gzip-compressed. The check after a tool that exited 0: it catches
# a run that wrote nothing where the tool said it succeeded.
wrote_vcf() {
  local first
  [ -s "$1" ] || return 1
  first=$(gzip -cdf "$1" 2>/dev/null | head -c 16) || true
  [[ "$first" == "##fileformat=VCF"* ]]
}

# vcf_header_samples: the sample names of the VCF header on stdin, one per
# line: the #CHROM columns after FORMAT. Nothing for a sites-only VCF. The
# pipeline analyses one sample per run (VCF_PRECHECK applies the same rule).
vcf_header_samples() {
  awk -F'\t' '/^#CHROM/ { for (i = 10; i <= NF; i++) print $i }'
}

# _bgzf_complete FILE: true when FILE ends with the 28-byte BGZF EOF block.
_bgzf_complete() {
  local last
  last=$(tail -c 28 "$1" | od -An -v -tx1 | tr -d ' \n')
  [ "$last" = "1f8b08040000000000ff0600424302001b0003000000000000000000" ]
}

# atomic_out FILE COMMAND [ARGS...]: run COMMAND with its standard output in
# FILE.tmp, and rename that to FILE only when COMMAND exits 0 and wrote
# something. Otherwise FILE.tmp is removed, FILE is left as it was and the
# exit code is returned, so a failed image pull or a crash never leaves an
# empty FILE behind for the next run to accept.
atomic_out() {
  local out=$1 rc=0
  shift
  rm -f "${out}.tmp"
  "$@" > "${out}.tmp" || rc=$?
  if [ "$rc" -eq 0 ] && [ ! -s "${out}.tmp" ]; then
    echo "ERROR: $1 wrote nothing for ${out}." >&2
    rc=1
  fi
  if [ "$rc" -ne 0 ]; then
    rm -f "${out}.tmp"
    return "$rc"
  fi
  mv -f "${out}.tmp" "$out"
}

# _digest md5|sha256|sum FILE: the checksum of FILE (sum = BSD sum, the
# format of Ensembl's CHECKSUMS files; only its first field is printed).
_digest() {
  case "$1" in
    md5)
      if command -v md5sum >/dev/null 2>&1; then md5sum "$2" | awk '{print $1}'
      elif command -v md5 >/dev/null 2>&1; then md5 -q "$2"
      else openssl md5 -r "$2" | awk '{print $1}'; fi ;;
    sha256)
      if command -v sha256sum >/dev/null 2>&1; then sha256sum "$2" | awk '{print $1}'; else shasum -a 256 "$2" | awk '{print $1}'; fi ;;
    sum) sum "$2" | awk '{print $1 + 0}' ;;
    *) echo "ERROR: unknown checksum type '${1}'" >&2; return 2 ;;
  esac
}

# consumer_sites_digest TABLE: the sha256 of the step-40 site table's content,
# read through gzip. The gzip wrapper holds the time it was written, so two
# builds of the same table differ as files and agree here. setup.sh and
# validate-setup.sh compare it with CONSUMER_SITES_TABLE_SHA256.
consumer_sites_digest() {
  _digest sha256 <(gzip -dc "$1")
}

# _get URL OUT: download URL to OUT (resuming a partial OUT), curl first and
# wget when curl is missing. OUT "-" writes to stdout. A server that cannot
# resume (curl exit 33) gets the partial OUT dropped and one fresh download.
#
# A stalled transfer fails instead of hanging: curl gives up on a connection
# that takes over 30 s, or a transfer slower than 1 KB/s for FETCH_STALL
# seconds (default 120); wget on a connection over 30 s, or no data for
# FETCH_STALL seconds, with no retry of its own. fetch then tries again,
# resuming the partial OUT.
_get() {
  local rc=0 stall=${FETCH_STALL:-120}
  # 10# so 00, 08 and 09 read as decimal; curl and wget take 0 as no bound.
  case "$stall" in ''|*[!0-9]*) stall=0 ;; esac
  [ "$((10#$stall))" -gt 0 ] || { echo "ERROR: FETCH_STALL must be a whole number of seconds above 0, got '${FETCH_STALL}'" >&2; return 1; }
  # The URL stays the last argument (scripts/ci/fake-docker/curl takes it from there).
  local cbound=(--connect-timeout 30 --speed-limit 1024 --speed-time "$stall")
  if command -v curl >/dev/null 2>&1; then
    if [ "$2" = "-" ]; then
      curl -fsSL "${cbound[@]}" "$1"
    else
      curl -fL "${cbound[@]}" -C - -o "$2" "$1" || rc=$?
      if [ "$rc" -eq 33 ]; then
        rm -f "$2"
        curl -fL "${cbound[@]}" -o "$2" "$1"
      else
        return "$rc"
      fi
    fi
  elif command -v wget >/dev/null 2>&1; then
    # --tries=1: fetch does the retries, as it does for curl.
    local wbound=(--connect-timeout=30 "--read-timeout=${stall}" --tries=1)
    if [ "$2" = "-" ]; then wget -q "${wbound[@]}" -O - "$1"; else wget -c "${wbound[@]}" -O "$2" "$1"; fi
  else
    echo "ERROR: neither curl nor wget is installed." >&2
    return 1
  fi
}

# fetch URL DEST [md5|sha256|sum VALUE|CHECKSUM_URL]
#
# Downloads URL to DEST.part (resuming one left by an interrupted run), checks
# it, then renames it to DEST. Nothing ever writes DEST directly, so a DEST
# that exists is whole. The check is the checksum when one is given (VALUE, or
# the file at CHECKSUM_URL: its only line, or the line naming URL's file), else
# `gzip -t` for a .gz file, else non-empty.
#
# Each try downloads, then checks. A transfer that broke keeps DEST.part, and
# the next try resumes it. A download that fails its check is deleted, and the
# next try downloads it from the start: one bad transfer, or a bad .part left
# by an earlier run, does not fail the call. The file and its CHECKSUM_URL are
# each tried FETCH_TRIES times (default 3), FETCH_WAIT seconds apart (default
# 5, then 10, ...); after the last try fetch returns 1. A larger FETCH_TRIES
# and FETCH_WAIT wait out a server that answers 404 for minutes at a time and
# then comes back.
fetch() {
  local url=$1 dest=$2 kind=${3:-} want=${4:-} part="${2}.part" name line i rc bad=""
  local tries=${FETCH_TRIES:-3} wait=${FETCH_WAIT:-}
  case "$tries" in ''|*[!0-9]*) tries=0 ;; *) tries=$((10#$tries)) ;; esac
  [ "$tries" -gt 0 ] || { echo "ERROR: FETCH_TRIES must be a whole number above 0, got '${FETCH_TRIES}'" >&2; return 1; }
  case "$wait" in *[!0-9]*) echo "ERROR: FETCH_WAIT must be a whole number of seconds, got '${wait}'" >&2; return 1 ;; esac
  name=$(basename "$url")
  if [ -n "$kind" ]; then
    case "$want" in
      http://*|https://*|ftp://*|file://*)
        # A read that worked but names no checksum is not retried.
        for ((i = 1; i <= tries; i++)); do
          if line=$(_get "$want" - | awk -v n="$name" '
                NF { lines++; first = $1; f = $NF; sub(/.*\//, "", f); if (f == n) hit = $1 }
                END { if (lines == 1) print first; else if (hit != "") print hit }'); then
            break
          fi
          line=""
          echo "  Checksum download attempt ${i}/${tries} failed: ${want}" >&2
          [ "$i" -lt "$tries" ] && sleep "${wait:-$((i * 5))}"
        done
        if [ -z "$line" ]; then
          echo "ERROR: could not read the checksum of ${name} from ${want}" >&2
          return 1
        fi
        want=$line ;;
    esac
  fi
  mkdir -p "$(dirname "$dest")"
  for ((i = 1; i <= tries; i++)); do
    if _get "$url" "$part"; then
      rc=0
      bad=$(_fetch_check "$kind" "$want" "$part" "$dest" "$name") || rc=$?
      [ "$rc" -eq 0 ] && { mv -f "$part" "$dest"; return; }
      [ "$rc" -eq 2 ] && return 1
      rm -f "$part"
      echo "  Download attempt ${i}/${tries} failed its check: ${bad}. Removed it." >&2
    else
      bad=""
      echo "  Download attempt ${i}/${tries} failed: ${url}" >&2
    fi
    [ "$i" -lt "$tries" ] && sleep "${wait:-$((i * 5))}"
  done
  if [ -n "$bad" ]; then
    echo "ERROR: ${bad}. Removed the download after ${tries} tries; run again to fetch it anew." >&2
  else
    echo "ERROR: could not download ${url} (the partial file ${part} is kept and resumed next time)" >&2
  fi
  return 1
}

# _fetch_check KIND WANT PART DEST NAME: fetch's check of one download. Exit 0
# when PART is good; 1 with the reason on stdout when it is not, and fetch
# downloads it again; 2 for an error a new download cannot fix.
_fetch_check() {
  local kind=$1 want=$2 part=$3 dest=$4 name=$5 got
  if [ -n "$kind" ]; then
    got=$(_digest "$kind" "$part") || return 2
    if [ "$kind" = sum ]; then
      case "$want" in
        ''|*[!0-9]*) echo "ERROR: '${want}' is not a sum checksum for ${name}" >&2; return 2 ;;
      esac
      want=$((10#${want}))
    fi
    if [ "$got" != "$want" ]; then
      echo "${kind} checksum of ${name} is ${got}, expected ${want}"
      return 1
    fi
  elif [[ "$dest" == *.gz ]]; then
    if ! gzip -t "$part" 2>/dev/null; then
      echo "${name} is not a complete gzip file"
      return 1
    fi
  elif [ ! -s "$part" ]; then
    echo "${name} downloaded empty"
    return 1
  fi
}

# --- Pinned data files ----------------------------------------------------------
# Small reference files from a fixed commit or release, each checked before it
# is stored. setup.sh installs them (install_data_file); a step only reads
# them (data_file) and runs no download of its own. HLA_DB_RELEASE can be
# overridden from the environment; step 08 keys its T1K index on it, so a new
# release builds a new index.
HLA_DB_RELEASE=${HLA_DB_RELEASE:-3.66.0}
GENCODE_RELEASE=50
DELLY_EXCL_COMMIT=72a05f3cd23bf55f1edd723b9504a0faf1fba230
# shellcheck disable=SC2034  # read by setup.sh and validate-setup.sh
DATA_FILES="delly_exclude cytoband hla_dat gencode_genes gencode_exons"

# data_file NAME: print the host path of data file NAME and return 0 when it
# is installed, 1 when it is not (the path is printed either way).
#   delly_exclude  Delly's GRCh38 exclude map: telomeres, centromeres and the
#                  extra contigs (step 19)
#   cytoband       UCSC's GRCh38 chromosome bands, chr1-22, X and Y (step 10)
#   hla_dat        hla.dat of IPD-IMGT/HLA release HLA_DB_RELEASE (step 08)
#   gencode_genes  the gene lines of GENCODE's basic annotation, release
#                  GENCODE_RELEASE: T1K takes HLA gene coordinates from it (step 08)
#   gencode_exons  the exons of the same annotation as a BED (chrom, 0-based
#                  start, end, once each): the cDNA read check of steps 06 and 17
data_file() {
  local out
  case "$1" in
    delly_exclude) out="${GENOME_DIR:?}/reference/delly_human.hg38.excl.tsv" ;;
    cytoband) out="${GENOME_DIR:?}/reference/cytoBand.hg38.txt" ;;
    hla_dat) out="${GENOME_DIR:?}/hla/IPD-IMGT-HLA_${HLA_DB_RELEASE}/hla.dat" ;;
    gencode_genes) out="${GENOME_DIR:?}/reference/gencode.v${GENCODE_RELEASE}.basic.genes.gtf" ;;
    gencode_exons) out="${GENOME_DIR:?}/reference/gencode.v${GENCODE_RELEASE}.basic.exons.bed" ;;
    *) echo "ERROR: unknown data file '${1}'" >&2; return 2 ;;
  esac
  printf '%s\n' "$out"
  [ -s "$out" ]
}

# install_data_file NAME: download, check and store data file NAME unless it
# is installed. Progress goes to stderr. Returns 1 when it cannot.
install_data_file() {
  local out url gz md5
  out=$(data_file "$1") && return 0
  [ -n "$out" ] || return 2
  case "$1" in
    delly_exclude)
      fetch "https://raw.githubusercontent.com/dellytools/delly/${DELLY_EXCL_COMMIT}/excludeTemplates/human.hg38.excl.tsv" \
        "$out" sha256 caef8593f82f513694ae41b429680eab8abe39c057a119dae7760715aef006ee >&2 || return 1 ;;
    cytoband)
      gz="${out}.download.gz"
      fetch https://hgdownload.soe.ucsc.edu/goldenPath/hg38/database/cytoBand.txt.gz "$gz" \
        sha256 e514e48b0a12a516fd5231b44241a20a282f290bdd6131cadcddc6276ee26082 >&2 || return 1
      # TelomereHunter reads chr1-22, X and Y and reports every other line
      # (the ALT, random and chrM rows) as invalid.
      _keep_lines "$gz" "$out" '$1 ~ /^chr([0-9]+|X|Y)$/' || return 1 ;;
    hla_dat)
      url="https://raw.githubusercontent.com/ANHIG/IMGTHLA/v${HLA_DB_RELEASE}-alpha"
      # md5checksum.txt lines read "MD5 (hla.dat.zip) = <md5>".
      md5=$(_get "${url}/md5checksum.txt" - 2>/dev/null | awk '$2 == "(hla.dat.zip)" {print $NF}') || md5=""
      if [ -z "$md5" ]; then
        echo "ERROR: no md5 for hla.dat.zip in ${url}/md5checksum.txt (is ${HLA_DB_RELEASE} an IPD-IMGT/HLA release?)" >&2
        return 1
      fi
      fetch "${url}/hla.dat.zip" "${out}.zip" md5 "$md5" >&2 || return 1
      # The T1K image has Perl's unzip module; the host may have no unzip.
      # shellcheck disable=SC2154  # T1K_IMAGE comes from versions.env
      run_in --rw "$(dirname "$out")" "$T1K_IMAGE" perl -e \
        'use IO::Uncompress::Unzip qw(unzip $UnzipError); unzip($ARGV[0] => $ARGV[1]) or die "unzip failed: $UnzipError\n";' \
        "$(cpath "${out}.zip")" "$(cpath "${out}.tmp")" >&2 || { rm -f "${out}.tmp"; return 1; }
      if ! grep -q 'IPD-IMGT/HLA' "${out}.tmp" 2>/dev/null; then
        rm -f "${out}.tmp"
        echo "ERROR: ${out}.zip held no IPD-IMGT/HLA hla.dat" >&2
        return 1
      fi
      mv -f "${out}.tmp" "$out"
      rm -f "${out}.zip" ;;
    gencode_genes|gencode_exons)
      _install_gencode || return 1 ;;
  esac
  [ -s "$out" ]
}

# _install_gencode: download GENCODE's basic annotation once and write from it
# each of gencode_genes and gencode_exons that is missing, then delete the
# download. Fails when one of them cannot be written.
_install_gencode() {
  local genes exons url gz rc=0
  genes=$(data_file gencode_genes) || true
  exons=$(data_file gencode_exons) || true
  url="https://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_human/release_${GENCODE_RELEASE}/gencode.v${GENCODE_RELEASE}.basic.annotation.gtf.gz"
  gz="${genes}.download.gz"
  fetch "$url" "$gz" md5 "$(dirname "$url")/MD5SUMS" >&2 || return 1
  if [ ! -s "$genes" ]; then
    _gz_to "$gz" "$genes" awk -F'\t' '$3 == "gene"' || rc=1
  fi
  if [ ! -s "$exons" ]; then
    _gz_to "$gz" "$exons" _gtf_exons || rc=1
  fi
  rm -f "$gz"
  return "$rc"
}

# _gtf_exons: the exon lines of a GTF on stdin as a sorted BED, once each.
_gtf_exons() {
  awk -F'\t' -v OFS='\t' '$3 == "exon" {print $1, $4 - 1, $5}' | LC_ALL=C sort -k1,1 -k2,2n -k3,3n -u
}

# _gz_to GZ OUT COMMAND...: write GZ, uncompressed and piped through COMMAND,
# to OUT (through OUT.tmp). Fails, keeping no OUT, when COMMAND writes nothing.
_gz_to() {
  local gz=$1 out=$2
  shift 2
  if gzip -cd "$gz" | "$@" > "${out}.tmp" && [ -s "${out}.tmp" ]; then
    mv -f "${out}.tmp" "$out"
  else
    rm -f "${out}.tmp"
    echo "ERROR: no usable lines for $(basename "$out") in $(basename "$gz")" >&2
    return 1
  fi
}

# _keep_lines GZ OUT AWK_CONDITION: write the tab-separated lines of GZ that
# match AWK_CONDITION to OUT (through OUT.tmp) and remove GZ. Fails, keeping
# no OUT, when no line matches.
_keep_lines() {
  if gzip -cd "$1" | awk -F'\t' "$3" > "${2}.tmp" && [ -s "${2}.tmp" ]; then
    mv -f "${2}.tmp" "$2"
    rm -f "$1"
  else
    rm -f "${2}.tmp" "$1"
    echo "ERROR: no usable lines in $(basename "$1")" >&2
    return 1
  fi
}

# vep_cache_url RELEASE: where Ensembl publishes the indexed GRCh38 VEP cache.
vep_cache_url() {
  printf 'https://ftp.ensembl.org/pub/release-%s/variation/indexed_vep_cache/homo_sapiens_vep_%s_GRCh38.tar.gz' "$1" "$1"
}

# install_vep_cache DIR RELEASE: download the cache (checked against Ensembl's
# CHECKSUMS file), unpack it in a temporary directory and move it to
# DIR/homo_sapiens/RELEASE_GRCh38 only when complete. The tarball is deleted.
install_vep_cache() {
  local dir=$1 rel=$2 url tarball tmp final
  url=$(vep_cache_url "$rel")
  tarball="${dir}/$(basename "$url")"
  final="${dir}/homo_sapiens/${rel}_GRCh38"
  if [ -e "$final" ]; then
    echo "ERROR: ${final} exists but has no info.txt (an interrupted extraction?)." >&2
    echo "  Remove it and run again." >&2
    return 1
  fi
  mkdir -p "${dir}/homo_sapiens"
  fetch "$url" "$tarball" sum "$(dirname "$url")/CHECKSUMS" || return 1
  echo "  Extracting $(basename "$tarball")..."
  tmp=$(mktemp -d "${dir}/.extract.XXXXXX")
  if ! tar xzf "$tarball" -C "$tmp" || [ ! -f "${tmp}/homo_sapiens/${rel}_GRCh38/info.txt" ]; then
    rm -rf "$tmp"
    echo "ERROR: could not unpack ${tarball} into a ${rel}_GRCh38 cache." >&2
    return 1
  fi
  mv "${tmp}/homo_sapiens/${rel}_GRCh38" "$final"
  rm -rf "$tmp" "$tarball"
}

# lock_acquire DIR / lock_release DIR: a lock shared by concurrent runs. The
# lock is a directory holding the owner's PID. It is stale when that PID is
# gone (a killed run), or when it has no PID a minute after it was made (a run
# killed between mkdir and writing the PID). A stale lock is renamed before it
# is removed, so of two waiters that judged it stale only one takes it over.
lock_acquire() {
  local lock=$1 pid stale
  until mkdir "$lock" 2>/dev/null; do
    pid=$(cat "${lock}/pid" 2>/dev/null || true)
    stale=false
    if [ -n "$pid" ]; then
      kill -0 "$pid" 2>/dev/null || stale=true
    elif [ -n "$(find "$lock" -maxdepth 0 -mmin +1 2>/dev/null)" ]; then
      stale=true
    fi
    if $stale; then
      if mv "$lock" "${lock}.stale.$$" 2>/dev/null; then
        if [ -n "$pid" ]; then
          echo "  Removing stale lock ${lock} (process ${pid} is gone)."
        else
          echo "  Removing stale lock ${lock} (no owner after a minute)."
        fi
        rm -rf "${lock}.stale.$$"
      fi
      continue
    fi
    echo "  Waiting for ${lock} (held by process ${pid:-starting})..."
    sleep 5
  done
  echo "$$" > "${lock}/pid"
}
lock_release() {
  rm -rf "$1"
}

# pipeline_images: the value of every *_IMAGE in versions.env, one per line,
# except those whose line is marked `# optional` (pulled on first use).
pipeline_images() {
  local line name
  while IFS= read -r line || [ -n "$line" ]; do
    case "$line" in
      *"# optional"*) continue ;;
    esac
    if [[ "$line" =~ ^([A-Z0-9_]+_IMAGE)= ]]; then
      name=${BASH_REMATCH[1]}
      printf '%s\n' "${!name}"
    fi
  done < "${PGP_ROOT}/versions.env"
}

# pinned_tag IMAGE: the local tag setup.sh gives an image versions.env pins by
# digest (name[:tag]@sha256:HEX), printed as localhost/pgp-pinned/NAME:sha256-HEX;
# nothing for an image pinned by tag alone. An image pulled by digest gets no
# tag, so Docker's classic image store lists it as dangling, and
# `docker image prune` or `docker rmi $(docker images -qf dangling=true)`
# deletes it. This second name keeps it out of that list. (The containerd
# image store keeps the digest reference as the image's name, so it never lists
# such an image as dangling.) The tag does not stop `docker system prune -a` or
# an explicit `docker rmi`. Every task still runs the digest reference: the tag
# only names the same image id. `localhost` is a registry host no one else can
# own, so the name can never be pulled from, or pushed to, a public registry.
pinned_tag() {
  local name hex last
  [[ "$1" =~ ^([a-z0-9][^@]*)@sha256:([0-9a-f]{64})$ ]] || return 0
  name=${BASH_REMATCH[1]} hex=${BASH_REMATCH[2]}
  last=${name##*/}
  [[ "$last" == *:* ]] && name=${name%:*}   # drop the tag, keep a registry port
  printf 'localhost/pgp-pinned/%s:sha256-%s\n' "${name//:/-}" "$hex"
}

# tag_pinned_image IMAGE: give a pulled digest-pinned image its pinned_tag.
# Safe to repeat: tagging an image with the name it already has changes nothing.
tag_pinned_image() {
  local tag
  tag=$(pinned_tag "$1")
  [ -z "$tag" ] || "$CONTAINER_ENGINE" tag "$1" "$tag"
}

# pinned_tag_present IMAGE: true when IMAGE needs no pinned_tag (it is pinned
# by tag alone) or its pinned_tag already names the same image id.
pinned_tag_present() {
  local tag
  tag=$(pinned_tag "$1")
  [ -z "$tag" ] ||
    [ "$("$CONTAINER_ENGINE" image inspect --format '{{.Id}}' "$tag" 2>/dev/null)" = \
      "$("$CONTAINER_ENGINE" image inspect --format '{{.Id}}' "$1" 2>/dev/null)" ]
}

# image_is_dangling IMAGE: 0 when Docker lists IMAGE's id as dangling (no tag
# of any name), the images `docker image prune` deletes; 1 when it does not
# (another local tag on the same id, the pinned_tag or not, or the containerd
# image store, where the digest reference is the image's name); 2 when
# Docker gives no id or no dangling list, so the answer is unknown.
image_is_dangling() {
  local id dangling
  id=$("$CONTAINER_ENGINE" image inspect --format '{{.Id}}' "$1" 2>/dev/null) && [ -n "$id" ] || return 2
  dangling=$("$CONTAINER_ENGINE" images --no-trunc --quiet --filter dangling=true 2>/dev/null) || return 2
  grep -qxF "$id" <<< "$dangling"
}

# scatter_beds DIR [INTERVALS]: split the calling of a single-process caller
# into units, one BED file each (DIR/001.bed, ...), and print their paths in
# reference order. Without INTERVALS: chr1-22, chrX, chrY and chrM of
# REF_FASTA's .fai one unit each, and every other contig in one last unit.
# With INTERVALS (space-separated contigs or contig:start-end regions, 1-based
# inclusive): one unit per region. SCATTER=false puts every region in one unit,
# so the caller runs once over all of them.
scatter_beds() {
  local dir=$1 intervals=${2:-} fai="${REF_FASTA}.fai" r
  rm -rf "$dir"
  mkdir -p "$dir"
  [ -s "$fai" ] || { echo "ERROR: ${fai} not found" >&2; return 1; }
  {
    if [ -n "$intervals" ]; then
      for r in $intervals; do printf 'R\t%s\n' "$r"; done
    fi
  } | awk -F'\t' -v d="$dir" -v all="${SCATTER:-true}" -v given="${intervals:+1}" '
    FNR == NR { len[$1] = $2; order[++n] = $1; rank[$1] = n; next }
    function unit(name) { if (all == "false") return sprintf("%s/001.bed", d); return sprintf("%s/%03d.bed", d, ++u) }
    { r = $2; gsub(/,/, "", r); c = r; s = 0; e = ""
      if (match(r, /:[0-9]+-[0-9]+$/)) { c = substr(r, 1, RSTART - 1); split(substr(r, RSTART + 1), p, "-"); s = p[1] - 1; e = p[2] }
      if (!(c in len)) { print "ERROR: contig " c " (INTERVALS) is not in the reference" > "/dev/stderr"; bad = 1; exit 1 }
      if (e == "") e = len[c]
      k++; rk[k] = rank[c]; st[k] = s + 0; reg[k] = c "\t" s "\t" e }
    END {
      if (bad) exit 1
      if (given) {
        # Reference order, then start, whatever the order of INTERVALS: the
        # units are joined in this order.
        for (i = 2; i <= k; i++)
          for (j = i; j > 1 && (rk[j - 1] > rk[j] || (rk[j - 1] == rk[j] && st[j - 1] > st[j])); j--) {
            t = rk[j]; rk[j] = rk[j - 1]; rk[j - 1] = t
            t = st[j]; st[j] = st[j - 1]; st[j - 1] = t
            t = reg[j]; reg[j] = reg[j - 1]; reg[j - 1] = t
          }
        for (i = 1; i <= k; i++) { f = unit(); print reg[i] >> f; close(f) }
        exit 0
      }
      for (i = 1; i <= n; i++) if (order[i] ~ /^chr([0-9]+|X|Y|M)$/) { f = unit(); print order[i] "\t0\t" len[order[i]] >> f; close(f) }
      rest = (all == "false") ? sprintf("%s/001.bed", d) : sprintf("%s/%03d.bed", d, u + 1)
      for (i = 1; i <= n; i++) if (order[i] !~ /^chr([0-9]+|X|Y|M)$/) print order[i] "\t0\t" len[order[i]] >> rest
    }' "$fai" - || return 1
  ls "$dir"/*.bed
}

# run_parallel MAX FUNC ARG...: run `FUNC ARG` for every ARG in the
# background, at most MAX at a time, and return 1 when any of them failed,
# after all have ended. Each unit's status comes from `wait PID`, never from
# `wait -n`: in CI (bash 5.2) `wait -n` sometimes returned 127, "no job to
# wait for", where a unit had ended, and that unit was counted as failed.
# Here `wait -n` only throttles and its status is not counted.
# Needs bash 4.3 (wait -n); the step scripts need 4.4.
run_parallel() {
  local max=$1 fn=$2 a p rc=0
  local pids=() live=() still=()
  shift 2
  for a in "$@"; do
    # Start the next unit once fewer than MAX are alive. `live` drops a PID
    # the first time `kill -0` fails, so a PID the kernel later reuses for
    # another process never holds a slot. A 127 from `wait -n` ends the wait
    # rather than looping on it.
    while :; do
      still=()
      for p in ${live[@]+"${live[@]}"}; do
        if kill -0 "$p" 2>/dev/null; then still+=("$p"); fi
      done
      live=(${still[@]+"${still[@]}"})
      [ "${#live[@]}" -ge "$max" ] || break
      wait -n || [ $? -ne 127 ] || break
    done
    "$fn" "$a" &
    pids+=("$!")
    live+=("$!")
  done
  for p in ${pids[@]+"${pids[@]}"}; do
    wait "$p" || rc=1
  done
  return "$rc"
}

# live_post RUN JSON: send one message to the live site (bin/owngenome_live.py,
# docs/live-site.md) at OWNGENOME_LIVE_URL/ingest/RUN; nothing without OWNGENOME_LIVE_URL.
# A site that is down or slow (3 s) costs one warning per script and never the
# step: it always returns 0.
live_post() {
  [ -n "${OWNGENOME_LIVE_URL:-}" ] || return 0
  if printf '%s' "$2" | curl -fsS -m 3 -o /dev/null -H 'Content-Type: application/json' -d @- \
      "${OWNGENOME_LIVE_URL%/}/ingest/$1" 2>/dev/null; then
    return 0
  fi
  if [ -z "${OWNGENOME_LIVE_WARNED:-}" ]; then
    echo "WARNING: the live site at ${OWNGENOME_LIVE_URL} did not take a message; the run goes on without it" >&2
    OWNGENOME_LIVE_WARNED=1
  fi
  return 0
}