From 89e0ae15a696ca8d8b1f2c5c2f28155a6504ca12 Mon Sep 17 00:00:00 2001 From: sphemakh Date: Wed, 11 Jun 2025 15:53:53 +0200 Subject: [PATCH 1/3] first steps --- .gitignore | 7 +++++++ contsub/parser/dimcontsub.py | 18 +++++++++++++++--- contsub/parser/imcontsub.yaml | 6 +++++- 3 files changed, 27 insertions(+), 4 deletions(-) diff --git a/.gitignore b/.gitignore index b1cb160..9b812a1 100644 --- a/.gitignore +++ b/.gitignore @@ -3,6 +3,13 @@ __pycache__/ *.py[cod] *$py.class +# Radio Astro +*.fits +*.zarr +*.ms +*.ipynb + + # C extensions *.so diff --git a/contsub/parser/dimcontsub.py b/contsub/parser/dimcontsub.py index d89a6c7..3788670 100644 --- a/contsub/parser/dimcontsub.py +++ b/contsub/parser/dimcontsub.py @@ -13,7 +13,6 @@ import dask.array as da import time import numpy as np -import xarray as xr import dask.multiprocessing log = init_logger(BIN.im_plane) @@ -64,11 +63,12 @@ def runit(**kwargs): else: cube = zds.DATA - sdim = (None,) # scalar dimension niter = 1 nomask = True + filemask = False if getattr(opts, "mask_image", None): mask = zds_from_fits(opts.mask_image, chunks=chunks).DATA + filemask = True nomask = False @@ -85,22 +85,34 @@ def runit(**kwargs): dask.config.set(scheduler='threads', num_workers = opts.nworkers) + prev_sclip = opts.sigma_clip[0] + sigma_clip = list(opts.sigma_clip) for iter_i in range(niter): futures = [] fitfunc = FitBSpline(opts.order[iter_i], opts.segments[iter_i]) + if nomask: + try: + sclip = sigma_clip[iter_i] + except IndexError: + sclip = prev_sclip + finally: + prev_sclip = sclip + for biter,dblock in enumerate(cube.data.blocks): if nomask and opts.automask: mask_future = get_mask(xspec, dblock, - opts.sigma_clip[0], + sclip, opts.order[iter_i], opts.segments[iter_i], ) + nomask = False elif nomask is False: mask_future = mask.data.blocks[biter] else: mask_future = da.zeros_like(dblock, dtype=bool) + contfit = ContSub(fitfunc, nomask=False, reshape=False, fitsaxes=False) getfit = da.gufunc( diff --git a/contsub/parser/imcontsub.yaml b/contsub/parser/imcontsub.yaml index 1dfdf0f..633ce01 100644 --- a/contsub/parser/imcontsub.yaml +++ b/contsub/parser/imcontsub.yaml @@ -27,7 +27,11 @@ inputs: default: - 5 - 4 - - 3 + - 3 + automask-per-iter: + info: Generate a new mask per iteration. + dtype: bool + default: False fit-model: dtype: str choices: From fa184aecedfb39ca86fc91ac3472f877900e24ec Mon Sep 17 00:00:00 2001 From: sphemakh Date: Fri, 29 Aug 2025 16:35:00 +0200 Subject: [PATCH 2/3] Fixes #10 --- contsub/parser/imcontsub.py | 8 ++++---- contsub/parser/imcontsub.yaml | 3 +++ contsub/utils.py | 4 +++- 3 files changed, 10 insertions(+), 5 deletions(-) diff --git a/contsub/parser/imcontsub.py b/contsub/parser/imcontsub.py index 4763c70..3c9d375 100644 --- a/contsub/parser/imcontsub.py +++ b/contsub/parser/imcontsub.py @@ -55,8 +55,8 @@ def runit(**kwargs): chunks = dict(ra = opts.ra_chunks or 64, dec=None, spectral=None) - - zds = zds_from_fits(infits, chunks=chunks) + rest_freq = opts.rest_freq + zds = zds_from_fits(infits, chunks=chunks, rest_freq=rest_freq) base_dims = ["ra", "dec", "spectral", "stokes"] if not hasattr(zds, "stokes"): base_dims.remove("stokes") @@ -75,13 +75,13 @@ def runit(**kwargs): if len(opts.order) != len(opts.segments): raise ValueError("If setting multiple --order and --segments, they must be of equal length. " - f"Got {len(opts.order)} orders and {len(opts.segments)} segments.") + f"Got {len(opts.order)} orders and {len(opts.segments)} segments.") niter = len(opts.order) nomask = True automask = False if getattr(opts, "mask_image", None): - mask = zds_from_fits(opts.mask_image, chunks=chunks).DATA + mask = zds_from_fits(opts.mask_image, chunks=chunks, rest_freq=rest_freq).DATA nomask = False if getattr(opts, "sigma_clip", None): diff --git a/contsub/parser/imcontsub.yaml b/contsub/parser/imcontsub.yaml index 5af944e..e79df72 100644 --- a/contsub/parser/imcontsub.yaml +++ b/contsub/parser/imcontsub.yaml @@ -53,6 +53,9 @@ inputs: dtype: int default: 0 choices: [0,1,2,3] + rest-freq: + info: Cube rest frequency in MHz. Will ignore the one in the FITS header if it exists. + dtype: float stokes-axis: info: "### DEPRECATED #### Set this flag if the input image has a stokes dimension. (Default is True)." default: yes diff --git a/contsub/utils.py b/contsub/utils.py index 3087548..c5c518b 100644 --- a/contsub/utils.py +++ b/contsub/utils.py @@ -43,7 +43,7 @@ def get_automask(xspec, cube, sigma_clip=5, order=3, segments=400): return mask -def zds_from_fits(fname, chunks=None): +def zds_from_fits(fname, chunks=None, rest_freq=None): """ Creates Zarr store from a FITS file. The resulting array has dimensions = RA, DEC, SPECTRAL[, STOKES] @@ -62,6 +62,8 @@ def zds_from_fits(fname, chunks=None): fds = xds_from_fits(fname)[0] header = fds.hdu.header + if rest_freq: + header["RESTFRQ"] = rest_freq * 1e6 # set it to Hz wcs = WCS(header, naxis="spectral stokes".split()) From 63c72a7639038f3e2c4c46ae14f95286a5f951f3 Mon Sep 17 00:00:00 2001 From: sphemakh Date: Fri, 29 Aug 2025 16:35:34 +0200 Subject: [PATCH 3/3] bump to 1.0.3 --- pyproject.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pyproject.toml b/pyproject.toml index 1cadac4..7abd4b3 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,6 +1,6 @@ [tool.poetry] name = "contsub" -version = "1.0.2" +version = "1.0.3" description = "Radio astronomy data continuum subtraction tools" authors = ["Amir Kazemi-Moridani, Sphesihle Makhathini, Mika Naidoo"] license = "MIT"