Skip to content

Commit

Permalink
prep and cal handles picket fence multichannel
Browse files Browse the repository at this point in the history
  • Loading branch information
d3v-null committed Sep 29, 2024
1 parent 24a4eae commit 69ce23b
Show file tree
Hide file tree
Showing 2 changed files with 149 additions and 107 deletions.
62 changes: 39 additions & 23 deletions demo/05_prep.sh
Original file line number Diff line number Diff line change
Expand Up @@ -16,9 +16,9 @@ export obsid=${obsid:-1341914000}
# RAW #
# ### #
# check for raw files
export raw_glob=${outdir}/${obsid}/raw/${obsid}_2\*.fits
if ! eval ls -1 $raw_glob >/dev/null; then
echo "raw not present: $raw_glob , try ${SCRIPT_BASE}/02_download.sh"
export raw_pattern=${outdir}/${obsid}/raw/${obsid}_2\*.fits
if ! eval ls -1 $raw_pattern >/dev/null; then
echo "raw not present: $raw_pattern , try ${SCRIPT_BASE}/02_download.sh"
exit 1
fi
# check for metafits files
Expand All @@ -41,22 +41,32 @@ export timeres_s=8 # time resolution to average to in seconds
export edgewidth_khz=80 # edge width to flag on each coarse channel in kHz

mkdir -p "${outdir}/${obsid}/prep"

# Prep_uvfits is where we tell birli to write our preprocessed visibility files to.
# However, If the channels in the observation are non-contiguous, birli will output
# multiple files add a channel suffix to the output file names.
# So to test if the uvfits files exist, we need to look for a pattern.
export prep_uvfits="${outdir}/${obsid}/prep/birli_${obsid}.uvfits"
export prepqa="${prep_uvfits%%.uvfits}_qa.json"
export prep_uvfits_pattern=${outdir}/${obsid}/prep/birli_${obsid}\*.uvfits

set -eu
if [[ ! -f $prep_uvfits ]]; then
# --provided-chan-ranges
birli ${birli_args:-} \
# since we don't expect the uvfits files to exist the first time around, 2>/dev/null silences the warning
if ! eval ls -1 $prep_uvfits_pattern 2>/dev/null; then
echo "running birli on $raw_pattern" \
$([[ -n "${edgewidth_khz:-}" ]] && echo " edge width ${edgewidth_khz}kHz") \
$([[ -n "${freqres_khz:-}" ]] && echo " freq res ${freqres_khz}kHz") \
$([[ -n "${timeres_s:-}" ]] && echo " time res ${timeres_s}s") \
;
/software/projects/mwaeor/dev/setonix/2024.05/software/linux-sles15-zen3/gcc-12.2.0/birli-main-soerdr3lmpn6rottzxzdxb2iehej73eu/bin/birli ${birli_args:-} \
-m "${metafits}" \
$([[ -n "${edgewidth_khz:-}" ]] && echo "--flag-edge-width ${edgewidth_khz}") \
$([[ -n "${freqres_khz:-}" ]] && echo "--avg-freq-res ${freqres_khz}") \
$([[ -n "${timeres_s:-}" ]] && echo "--avg-time-res ${timeres_s}") \
-u "${prep_uvfits}" \
$raw_glob \
$raw_pattern \
$@
# -M "${prep_uvfits%%.uvfits}.ms" \
else
echo "prep_uvfits $prep_uvfits exists, skipping birli"
echo "prep_uvfits $prep_uvfits_pattern exists, skipping birli"
fi

# ####### #
Expand All @@ -65,19 +75,25 @@ fi
# DEMO: use mwa_qa for quality analysis of preprocessed uvfits
# details: https://github.com/d3v-null/mwa_qa (my fork of https://github.com/Chuneeta/mwa_qa/ )

if [[ ! -f "$prepqa" ]]; then
run_prepvisqa.py $prep_uvfits $metafits --out $prepqa
else
echo "prepqa=$prepqa exists, skipping run_prepvisqa.py"
fi
# loop over all the preprocessed visibility files birli produced
eval ls -1 $prep_uvfits_pattern | while read -r prep_uvfits; do
export prep_uvfits
# obsid plus any channel suffix that birli might add.
export obs_ch=${prep_uvfits##*/birli_}
export obs_ch=${obs_ch%%.uvfits}

# store prepqa relative to this uvfits
export prepqa="${prep_uvfits%%.uvfits}_qa.json"
if [[ ! -f "$prepqa" ]]; then
run_prepvisqa.py $prep_uvfits $metafits --out $prepqa
fi

# DEMO: extract bad antennas from prepqa json with jq
# - both of the provided observations pass QA, so no bad antennas are reported
prep_bad_ants=$(jq -r $'.BAD_ANTS|join(" ")' $prepqa)
# DEMO: extract bad antennas from prepqa json with jq
prep_bad_ants=$(jq -r $'.BAD_ANTS|join(" ")' $prepqa)

# DEMO: plot the prep qa results
# - RMS plot: RMS of all autocorrelation values for each antenna
# - zscore:
plot_prepvisqa.py $prepqa --save --out ${prep_uvfits%%.uvfits}
# DEMO: plot the prep qa results
# - RMS plot: RMS of all autocorrelation values for each antenna
plot_prepvisqa.py $prepqa --save --out ${prep_uvfits%%.uvfits}

echo $obsid $prep_bad_ants
echo $obs_ch $prep_bad_ants
done
194 changes: 110 additions & 84 deletions demo/06_cal.sh
Original file line number Diff line number Diff line change
Expand Up @@ -23,31 +23,23 @@ if [[ ! -f "$metafits" ]]; then
fi

# check preprocessed visibility and qa files exist from previous steps
export prep_uvfits="${outdir}/${obsid}/prep/birli_${obsid}.uvfits"
if [ ! -f "$prep_uvfits" ]; then
echo "prep_uvfits=$prep_uvfits does not exist. trying 05_prep.sh"
# - birli adds a channel suffix when processing observations with non-contiguous coarse channels.
# - if the files we need are missing, then run 05_prep.
export prep_uvfits_pattern=${outdir}/${obsid}/prep/birli_${obsid}\*.uvfits
if ! eval ls -1 $prep_uvfits_pattern >/dev/null; then
echo "prep_uvfits $prep_uvfits_pattern does not exist. trying 05_prep.sh"
$SCRIPT_BASE/05_prep.sh
fi
export prepqa="${prep_uvfits%%.uvfits}_qa.json"
if [[ ! -f "$prepqa" ]]; then
echo "warning: prepqa=$prepqa does not exist. trying 05_prep.sh"
export prep_qa_pattern="${outdir}/${obsid}/prep/birli_${obsid}*_qa.json"
if ! eval ls -1 $prep_qa_pattern >/dev/null; then
echo "warning: prepqa $prep_qa_pattern does not exist. trying 05_prep.sh"
$SCRIPT_BASE/05_prep.sh
fi
if [[ -f "$prepqa" ]]; then
prep_bad_ants=$(jq -r '.BAD_ANTS|join(" ")' "$prepqa")
export prep_bad_ants
else
export prep_bad_ants=""
fi

# ### #
# CAL #
# ### #
# direction independent calibration with hyperdrive
mkdir -p "${outdir}/${obsid}/cal"
export hyp_soln="${outdir}/${obsid}/cal/hyp_soln_${obsid}.fits"
export cal_ms="${outdir}/${obsid}/cal/hyp_cal_${obsid}.ms"
export model_ms="${outdir}/${obsid}/cal/hyp_model_${obsid}.ms"

# DEMO: generate calibration solutions from preprocessed visibilities
# details: https://mwatelescope.github.io/mwa_hyperdrive/user/di_cal/intro.html
Expand All @@ -59,72 +51,106 @@ if [[ -n "${gpus:-}" ]]; then
dical_args=""
fi

mkdir -p "${outdir}/${obsid}/cal"
set -eu

if [[ ! -f "$hyp_soln" ]]; then
echo "calibrating with sourcelist $srclist"
hyperdrive di-calibrate ${dical_args:-} \
--data "$metafits" "$prep_uvfits" \
--source-list "$srclist" \
--outputs "$hyp_soln" \
$([[ -n "${prep_bad_ants:-}" ]] && echo --tile-flags $prep_bad_ants)
else
echo "hyp_soln $hyp_soln exists, skipping hyperdrive di-calibrate"
fi

# plot solutions file
if [[ ! -f "${hyp_soln%%.fits}_phases.png" ]]; then
hyperdrive solutions-plot \
-m "$metafits" \
--no-ref-tile \
--max-amp 1.5 \
--output-directory "${outdir}/${obsid}/cal" \
"$hyp_soln"
else
echo "phases_png ${hyp_soln%%.fits}_phases.png exists, skipping hyperdrive solutions-plot"
fi

# ###### #
# CAL QA #
# ###### #
# DEMO: use mwa_qa to check calibration solutions

export calqa="${hyp_soln%%.fits}_qa.json"

if [[ ! -f "$calqa" ]]; then
echo "running calqa on solutions $hyp_soln"
run_calqa.py --pol X --out "$calqa" "$hyp_soln" "$metafits"
else
echo "calqa $calqa exists, skipping run_calqa.py"
fi

# plot the cal qa results
plot_calqa.py "$calqa" --save --out "${hyp_soln%%.fits}"

# extract the percentage of channels that converged
cal_pct_nonconvg=$(jq -r $'.PERCENT_NONCONVERGED_CHS|tonumber|round' "$calqa")
export cal_pct_nonconvg

if [[ $cal_pct_nonconvg -ge 95 ]]; then
echo "calibration failed, $cal_pct_nonconvg% of channels did not converge. hint: try a different sky model in demo/00_env.sh"
exit 1
fi

# extract bad antennas from calqa json with jq
cal_bad_ants=$(jq -r $'.BAD_ANTS|join(" ")' "$calqa")
export cal_bad_ants

echo "deliberately disabling cal bad ants for the first round :)"
export cal_bad_ants=""

# apply calibration solutions to preprocessed visibilities
# details: https://mwatelescope.github.io/mwa_hyperdrive/user/solutions_apply/intro.html
if [[ ! -d "$cal_ms" ]]; then
hyperdrive apply ${apply_args:-} \
--data "$metafits" "$prep_uvfits" \
--solutions "$hyp_soln" \
--outputs "$cal_ms" \
$([[ -n "${cal_bad_ants:-}" ]] && echo --tile-flags $cal_bad_ants)
else
echo "cal_ms $cal_ms exists, skipping hyperdrive apply"
fi
# loop over all the preprocessed files
eval ls -1 $prep_uvfits_pattern | while read -r prep_uvfits; do
export prep_uvfits;

# find prepqa relative to this uvfits file
export prepqa="${prep_uvfits%%.uvfits}_qa.json"
if [[ -f "$prepqa" ]]; then
prep_bad_ants=$(jq -r '.BAD_ANTS|join(" ")' "$prepqa")
export prep_bad_ants
else
export prep_bad_ants=""
fi

# store calibration outputs for the prepqa file in the sibling cal/ folder
# e.g. for prep_uvfits=a/b/prep/birli_X_chY.uvfits, parent=a/b, obs=X_chY
export parent=${prep_uvfits%/*}
export parent=${parent%/*}
export obs=${prep_uvfits##*/birli_}
export obs=${obs%.uvfits}
export hyp_soln="${parent}/cal/hyp_soln_${obs}.fits"
export cal_ms="${parent}/cal/hyp_cal_${obs}.ms"
export model_ms="${parent}/cal/hyp_model_${obs}.ms"

if [[ ! -f "$hyp_soln" ]]; then
echo "calibrating with sourcelist $srclist"
hyperdrive di-calibrate ${dical_args:-} \
--data "$metafits" "$prep_uvfits" \
--source-list "$srclist" \
--outputs "$hyp_soln" \
$([[ -n "${prep_bad_ants:-}" ]] && echo --tile-flags $prep_bad_ants)
# TODO: --cpu if login node
else
echo "hyp_soln $hyp_soln exists, skipping hyperdrive di-calibrate"
fi

# plot solutions file
if [[ ! -f "${hyp_soln%%.fits}_phases.png" ]]; then
hyperdrive solutions-plot \
-m "$metafits" \
--no-ref-tile \
--max-amp 1.5 \
--output-directory "${outdir}/${obsid}/cal" \
"$hyp_soln"
else
echo "phases_png ${hyp_soln%%.fits}_phases.png exists, skipping hyperdrive solutions-plot"
fi

# ###### #
# CAL QA #
# ###### #
# DEMO: use mwa_qa to check calibration solutions

export calqa="${hyp_soln%%.fits}_qa.json"

if [[ ! -f "$calqa" ]]; then
echo "running calqa on solutions $hyp_soln"
run_calqa.py --pol X --out "$calqa" "$hyp_soln" "$metafits"
else
echo "calqa $calqa exists, skipping run_calqa.py"
fi

# plot the cal qa results
plot_calqa.py "$calqa" --save --out "${hyp_soln%%.fits}"

# extract the percentage of channels that converged
cal_pct_nonconvg=$(jq -r $'.PERCENT_NONCONVERGED_CHS|tonumber|round' "$calqa")
export cal_pct_nonconvg

if [[ $cal_pct_nonconvg -ge 95 ]]; then
echo "calibration failed, $cal_pct_nonconvg% of channels did not converge. hint: try a different sky model in demo/00_env.sh"
continue
fi

# extract bad antennas from calqa json with jq
cal_bad_ants=$(jq -r $'.BAD_ANTS|join(" ")' "$calqa")
export cal_bad_ants

# echo "deliberately disabling cal bad ants for the first round :)"
# export cal_bad_ants=""

# apply calibration solutions to preprocessed visibilities
# details: https://mwatelescope.github.io/mwa_hyperdrive/user/solutions_apply/intro.html
if [[ ! -d "$cal_ms" ]]; then
hyperdrive apply ${apply_args:-} \
--data "$metafits" "$prep_uvfits" \
--solutions "$hyp_soln" \
--outputs "$cal_ms" \
$([[ -n "${cal_bad_ants:-}" ]] && echo --tile-flags $cal_bad_ants)
else
echo "cal_ms $cal_ms exists, skipping hyperdrive apply"
fi
done

# TODO phase fits, for now you need to singularity exec -B$PWD -B${outdir:-$PWD} -W$PWD --cleanenv docker://mwatelescope/mwa-demo:latest /bin/bash
cat <<EOF
# this will only work in the contaier :(
python /mwax_mover/scripts/cal_analysis.py \
--name "${obsid}" \
--metafits "${metafits}" --solns ${outdir}/${obsid}/cal/hyp_soln_${obsid}*.fits \
--plot-residual --residual-vmax=0.5
EOF

0 comments on commit 69ce23b

Please sign in to comment.