From d1d33f706ea42dacf5cc113ebd0f3ed1f7fd79ef Mon Sep 17 00:00:00 2001 From: Thomas Williams <32936518+thomaswilliamsastro@users.noreply.github.com> Date: Wed, 16 Sep 2026 10:29:57 +0100 Subject: [PATCH] Use weight file from SD pipeline This PR brings in the weight files that the SD pipeline produces, to improve linear mosaicking. - Edit the weight file in the SD pipeline to square what the pipeline spits out, so the weight works as expected - Pull in the SD weight file, if it exists - Stage SD weight file like the SD file - Don't generate weight file if SD weight file already exists - Add explanation for .weight file in docs --- CHANGELOG.md | 1 + .../handlers/SingleDishHandler.rst | 5 +- phangsPipeline/casaSingleDishALMAWrapper.py | 9 +- phangsPipeline/handlerPostprocess.py | 118 +++++++++++------- 4 files changed, 87 insertions(+), 46 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index ea5d18c3..68b9ae5b 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -27,6 +27,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - Improved convolution with spectral-cube, including new analytic support for small kernels (#385). - Add tests for utilsLines, and tidy up formatting (#403). - Add optional file string in derived for nonstandard workflows (#406). +- Use SD pipeline weights in single-dish linear mosaicking (#397). ### Changed diff --git a/docs/docs_phangsPipeline/handlers/SingleDishHandler.rst b/docs/docs_phangsPipeline/handlers/SingleDishHandler.rst index f1a7a087..633ddd2a 100644 --- a/docs/docs_phangsPipeline/handlers/SingleDishHandler.rst +++ b/docs/docs_phangsPipeline/handlers/SingleDishHandler.rst @@ -4,7 +4,10 @@ SingleDishHandler The SingleDishHandler allows for imaging of singledish ALMA data. This essentially replicates what the ALMA pipeline does, and then -produces cubes at the end. +produces ``.image`` and ``.weight`` cubes at the end. The weight +cube scales roughly as :math:`t_int \times \Delta v / T_sys^2`, +and so are squared afterwards, so they can be fed directly +into later linear mosaicking. .. autofunction:: phangsPipeline.SingleDishHandler.loop_singledish :noindex: diff --git a/phangsPipeline/casaSingleDishALMAWrapper.py b/phangsPipeline/casaSingleDishALMAWrapper.py index d5ea31e4..50aee95e 100644 --- a/phangsPipeline/casaSingleDishALMAWrapper.py +++ b/phangsPipeline/casaSingleDishALMAWrapper.py @@ -6,6 +6,7 @@ import astropy.units as u import numpy as np +from astropy.io import fits from astropy.table import Table from . import casaLegacySingleDishRoutines as csdr @@ -407,10 +408,14 @@ def runALMAPipeline(path_galaxy, weightimage_fits = f"{weightimage}.fits" casaStuff.exportfits(imagename=weightimage, fitsimage=weightimage_fits, overwrite=True) + # Move the image shutil.copy2(imagefile_fits, output_file) - # And export the weightfile + + # Square the weight so this can be properly used in mosaicking later weight_output_file = output_file.replace(".fits", '_weight.fits') - shutil.copy2(weightimage_fits, weight_output_file) + with fits.open(weightimage_fits) as hdu: + hdu[0].data = hdu[0].data ** 2 + hdu.writeto(weight_output_file, overwrite=True) logger.info(f"> Copied image to {output_file}") logger.info(f"> Copied weight to {weight_output_file}") diff --git a/phangsPipeline/handlerPostprocess.py b/phangsPipeline/handlerPostprocess.py index 9b63b4b9..64463ba1 100644 --- a/phangsPipeline/handlerPostprocess.py +++ b/phangsPipeline/handlerPostprocess.py @@ -180,6 +180,16 @@ def _fname_dict( else: fname_dict[tag] = '' + tag = "orig_weight_sd" + if has_sd: + orig_weight_sd_file = self._kh.get_sd_filename( + target=target, + product=product, + ) + fname_dict[tag] = orig_weight_sd_file.replace(".fits", "_weight.fits") + else: + fname_dict[tag] = "" + # &%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&% # Processed files (apply the extra_ext tag here) # &%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&% @@ -657,59 +667,74 @@ def task_stage_singledish( target=target, config=config, product=product, extra_ext=extra_ext_out) template = fname_dict_in[template_tag] - infile = fname_dict_in['orig_sd'] - outfile = fname_dict_out[out_tag] - - # Check input file existence - if check_files: - - files_exist = check_files_exist(indir+infile) - if not files_exist: - return () - files_exist = check_files_exist(tempdir+template, postprocessing_method=postprocessing_method, ) if not files_exist: return () - # Stage the singledish data for feathering + for this_tag in [ + "orig_sd", + "orig_weight_sd", + ]: - logger.info("") - logger.info("&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%") - logger.info("Preparing single dish data for:") - logger.info(str(target) + " , " + str(product) + " , " + str(config)) - logger.info("&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%") - logger.info("") + out_tag = { + "orig_sd": "prepped_sd", + "orig_weight_sd": "sd_weight", + }[this_tag] - logger.info("Using prep_sd_for_feather.") - logger.info("Prepping " + outfile) - logger.info("Original file " + infile) - logger.info("Using interferometric template " + template) + infile = fname_dict_in[this_tag] + outfile = fname_dict_out[out_tag] - if not self._dry_run: - - if postprocessing_method == "casa": - cfr.prep_sd_for_feather( - sdfile_in=indir + infile, - sdfile_out=outdir + outfile, - interf_file=tempdir + template, - do_import=True, - do_align=True, - do_checkunits=True, - overwrite=True) - elif postprocessing_method == "spectralcube": - sfr.prep_sd_for_feather( - sdfile_in=f"{indir}{infile}", - sdfile_out=f"{outdir}{outfile}.fits", - interf_file=f"{tempdir}{template}.fits", - do_align=True, - do_checkunits=True, - overwrite=True, - ) - else: - raise ValueError(f"postprocessing_method must be one of {ALLOWED_POSTPROCESSING_METHODS}") + # Don't check units on the weight file + do_checkunits = True + if this_tag in ["orig_weight_sd"]: + do_checkunits = False + + # Check input file existence + if check_files: + + files_exist = check_files_exist(indir+infile) + if not files_exist: + return () + + # Stage the singledish data for feathering + + logger.info("") + logger.info("&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%") + logger.info("Preparing single dish data for:") + logger.info(str(target) + " , " + str(product) + " , " + str(config)) + logger.info("&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%&%") + logger.info("") + + logger.info("Using prep_sd_for_feather.") + logger.info("Prepping " + outfile) + logger.info("Original file " + infile) + logger.info("Using interferometric template " + template) + + if not self._dry_run: + + if postprocessing_method == "casa": + cfr.prep_sd_for_feather( + sdfile_in=indir + infile, + sdfile_out=outdir + outfile, + interf_file=tempdir + template, + do_import=True, + do_align=True, + do_checkunits=do_checkunits, + overwrite=True) + elif postprocessing_method == "spectralcube": + sfr.prep_sd_for_feather( + sdfile_in=f"{indir}{infile}", + sdfile_out=f"{outdir}{outfile}.fits", + interf_file=f"{tempdir}{template}.fits", + do_align=True, + do_checkunits=do_checkunits, + overwrite=True, + ) + else: + raise ValueError(f"postprocessing_method must be one of {ALLOWED_POSTPROCESSING_METHODS}") return () @@ -888,6 +913,13 @@ def task_make_singledish_weight( fname_dict_out = self._fname_dict( target=target, config=config, product=product, extra_ext=extra_ext_out) + # Check we haven't pulled in a singledish weight file already, and if we do, skip this step. + weight_file_exists = check_files_exist(indir+fname_dict_in["sd_weight"], + postprocessing_method=postprocessing_method, + ) + if weight_file_exists: + return () + image_file = fname_dict_in[image_tag] outfile = fname_dict_out[out_tag]