From 047928b7e2e2e20dfd1d08b6c41e56a753df89ab Mon Sep 17 00:00:00 2001 From: tgalvin Date: Thu, 29 Feb 2024 14:01:54 +0800 Subject: [PATCH 01/21] exposed preflagger args / mesh_ant_flags added --- flint/bptools/preflagger.py | 56 ++++++++++++++++++--- flint/calibrate/aocalibrate.py | 24 ++++++--- flint/options.py | 4 ++ flint/prefect/flows/bandpass_pipeline.py | 62 ++++++++++++------------ tests/test_aocalibrate.py | 21 ++++++++ 5 files changed, 125 insertions(+), 42 deletions(-) diff --git a/flint/bptools/preflagger.py b/flint/bptools/preflagger.py index 2d617314..39336c89 100644 --- a/flint/bptools/preflagger.py +++ b/flint/bptools/preflagger.py @@ -7,7 +7,7 @@ """ from pathlib import Path -from typing import NamedTuple, Optional, Tuple +from typing import NamedTuple, Optional, Tuple, List import matplotlib.pyplot as plt import numpy as np @@ -443,7 +443,7 @@ def flag_mean_residual_amplitude( def flag_mean_xxyy_amplitude_ratio( - xx_complex_gains: np.ndarray, yy_complex_gains, fraction: float = 2.0 + xx_complex_gains: np.ndarray, yy_complex_gains, tolerance: float = 0.1 ) -> bool: """Will robust compute through an iterative sigma-clipping procedure the mean XX and YY gain amplitudes. The ratio of these means are computed, @@ -457,7 +457,7 @@ def flag_mean_xxyy_amplitude_ratio( Args: xx_complex_gains (np.ndarray): The XX complex gains to be considered yy_complex_gains (_type_): The YY complex gains to be considered - fraction (float, optional): The fraction used to distinguish a critical mean ratio threshold. Defaults to 2.. + tolerance (float, optional): The tolerance used used to distinguish a critical mean ratio threshold. Defaults to 0.10. Returns: bool: Whether data should be flagged (True) or not (False) @@ -482,13 +482,57 @@ def flag_mean_xxyy_amplitude_ratio( result = ( not np.isfinite(mean_gain_ratio) - or mean_gain_ratio < (1.0 / fraction) - or mean_gain_ratio > fraction + or mean_gain_ratio < (1.0 - tolerance) + or mean_gain_ratio > (1.0 + tolerance) ) if result: logger.warning( - f"Failed the mean gain ratio test: {xx_mean=} {yy_mean=} {mean_gain_ratio=} " + f"Failed the mean gain ratio test: {xx_mean=} {yy_mean=} {mean_gain_ratio=} {tolerance=}" ) return result + + +def construct_mesh_ant_flags(mask: np.ndarray) -> np.ndarray: + """Construct a mask that will accumulate the flags across + all antennas. The input mask array should be boolean and + of shape (ant, channels, pol), where `True` means flagged. + + If an antenna is completely flagged it is ignored as the + statistics are collected + + Args: + mask (np.ndarray): Input array denoting which items are flagged. + + Returns: + np.ndarray: Output array where antennas have common sets of flags + """ + + assert ( + len(mask.shape) == 3 + ), f"Expect array of shape (ant, chnnel, pol), received {mask.shape=}" + accumulate_mask = np.zeros_like(mask[0], dtype=bool) + + nant = mask.shape[0] + logger.info(f"Accumulating flagged channels over {nant=} antenna") + + empty_ants: List[int] = [] + + for ant in range(nant): + ant_mask = mask[ant] + if np.all(ant_mask): + empty_ants.append(ant) + continue + + accumulate_mask = accumulate_mask | ant_mask + + result_mask = np.zeros_like(mask, dtype=bool) + + for ant in range(nant): + if ant in empty_ants: + result_mask[ant, :, :] = True + else: + result_mask[ant] = accumulate_mask + + return result_mask diff --git a/flint/calibrate/aocalibrate.py b/flint/calibrate/aocalibrate.py index 9e3d6f36..1ac63489 100644 --- a/flint/calibrate/aocalibrate.py +++ b/flint/calibrate/aocalibrate.py @@ -26,6 +26,7 @@ flag_mean_xxyy_amplitude_ratio, flag_outlier_phase, flags_over_threshold, + construct_mesh_ant_flags, ) from flint.bptools.smoother import ( divide_bandpass_by_ref_ant_preserve_phase, @@ -863,6 +864,8 @@ def flag_aosolutions( plot_solutions_throughout: bool = True, smooth_window_size: int = 16, smooth_polynomial_order: int = 4, + mean_ant_tolderance: float = 0.1, + mesh_ant_flags: bool = False, ) -> FlaggedAOSolution: """Will open a previously solved ao-calibrate solutions file and flag additional channels and antennae. @@ -884,6 +887,8 @@ def flag_aosolutions( plot_solutions_throughout (bool, Optional): If True, the solutions will be plotted at different stages of processing. Defaults to True. smooth_window_size (int, optional): The size of the window function of the savgol filter. Passed directly to savgol. Defaults to 16. smooth_polynomial_order (int, optional): The order of the polynomial of the savgol filter. Passed directly to savgol. Defaults to 4. + mean_ant_tolerance (float, optional): Tolerance of the mean x/y antenna gain ratio test before the antenna is flagged. Defaults to 0.1. + mesh_ant_flags (bool, optional): If True, a channel is flagged across all antenna if it is flagged for any antenna. Defaults to False. Returns: FlaggedAOSolution: Path to the updated solutions file, intermediate solution files and plots along the way @@ -919,6 +924,12 @@ def flag_aosolutions( output_plots = plot_solutions(solutions=solutions_path, ref_ant=ref_ant) plots.extend(output_plots) + if mesh_ant_flags: + logger.info("Combining antenna flags") + for time in range(solutions.nsol): + time_mask = construct_mesh_ant_flags(mask=~np.isfinite(bandpass[time])) + bandpass[time, time_mask] = np.nan + for time in range(solutions.nsol): for pol in (0, 3): logger.info(f"Processing {pols[pol]} polarisation") @@ -981,16 +992,17 @@ def flag_aosolutions( ) for time in range(solutions.nsol): - ref_ant_gains = bandpass[time, ref_ant] + bandpass_phased_referenced = divide_bandpass_by_ref_ant_preserve_phase( + complex_gains=bandpass[time], ref_ant=ref_ant + ) # This loop will flag based on stats across different polarisations for ant in range(solutions.nant): - # We need to skip the case of flagging on the reference antenna, I think. - if ref_ant == ant: - continue - ant_gains = bandpass[time, ant] / ref_ant_gains + ant_gains = bandpass_phased_referenced[ant] if flag_mean_xxyy_amplitude_ratio( - xx_complex_gains=ant_gains[:, 0], yy_complex_gains=ant_gains[:, 3] + xx_complex_gains=ant_gains[:, 0], + yy_complex_gains=ant_gains[:, 3], + tolerance=mean_ant_tolderance, ): logger.info(f"{ant=} failed mean amplitude gain test. Flagging {ant=}.") bandpass[time, ant, :, :] = np.nan diff --git a/flint/options.py b/flint/options.py index 2308297e..505e530b 100644 --- a/flint/options.py +++ b/flint/options.py @@ -33,6 +33,10 @@ class BandpassOptions(NamedTuple): """The number of times the bandpass will be calibrated, flagged, then recalibrated""" minuv: Optional[float] = None """The minimum baseline length, in meters, for data to be included in bandpass calibration stage""" + preflagger_ant_mean_tolerance: float = 0.2 + """Tolerance that the mean x/y antenna gain ratio test before the antenna is flagged""" + preflagger_mesh_ant_flags: bool = False + """Share channel flags from bandpass solutions between all antenna""" class FieldOptions(NamedTuple): diff --git a/flint/prefect/flows/bandpass_pipeline.py b/flint/prefect/flows/bandpass_pipeline.py index cf286ed9..30fb3d17 100644 --- a/flint/prefect/flows/bandpass_pipeline.py +++ b/flint/prefect/flows/bandpass_pipeline.py @@ -126,15 +126,10 @@ def task_flag_solutions( def run_bandpass_stage( bandpass_mss: Collection[MS], output_split_bandpass_path: Path, - calibrate_container: Path, - flagger_container: Path, + bandpass_options: BandpassOptions, model_path: Path, source_name_prefix: str = "B1934-638", skip_rotation: bool = False, - smooth_window_size: int = 16, - smooth_polynomial_order: int = 4, - flag_calibrate_rounds: int = 0, - minuv: Optional[float] = None, ) -> List[CalibrateCommand]: """Excutes the bandpass calibration (using ``calibrate``) against a set of input measurement sets. @@ -142,22 +137,17 @@ def run_bandpass_stage( Args: bandpass_mss (Collection[MS]): Set of bandpass measurement sets to calibrate output_split_bandpass_path (Path): The location where the extract field centred on the calibration field (typically PKSB19340638) - calibrate_container (Path): Path the a singularity container with the ``calibrate`` software - flagger_container (Path): Path to the singularity container with the ``aoflagger`` software + bandpass_options (BandpassOptions): Configurables that will specify the bandpass calibbration process model_path (Path): Path to the model used to calibrate against source_name_prefix (str, optional): Name of the field being calibrated. Defaults to "B1934-638". skip_rotation (bool, optional): If ``True`` the rotation of the ASKAP visibility from the antenna frame to the sky-frame will be skipped. Defaults to False. - smooth_window_size (int, optional): The size of the window function of the savgol filter. Passed directly to savgol. Defaults to 16. - smooth_polynomial_order (int, optional): The order of the polynomial of the savgol filter. Passed directly to savgol. Defaults to 4. - flag_calibrate_rounds (int, optional): Defines the length of a loop that will calibrate, apply, flag and recalibrate. If 0, this is not performed. Defaults to 0. - minuv (Optional[float], optional): The minimum baseline length, in meters, for data to be included in bandpass calibration stage. Defaults to None. Returns: List[CalibrateCommand]: Set of calibration commands used """ assert ( - flag_calibrate_rounds >= 0 - ), f"Currently {flag_calibrate_rounds=}, needs to be 0 or higher" + bandpass_options.flag_calibrate_rounds >= 0 + ), f"Currently {bandpass_options.flag_calibrate_rounds=}, needs to be 0 or higher" if not output_split_bandpass_path.exists(): logger.info(f"Creating {str(output_split_bandpass_path)}") @@ -174,37 +164,40 @@ def run_bandpass_stage( ms=extract_bandpass_mss, skip_rotation=skip_rotation ) flag_bandpass_mss = task_flag_ms_aoflagger.map( - ms=preprocess_bandpass_mss, container=flagger_container, rounds=1 + ms=preprocess_bandpass_mss, + container=bandpass_options.flagger_container, + rounds=1, ) calibrate_cmds = task_create_calibrate_cmd.map( ms=flag_bandpass_mss, calibrate_model=model_path, - container=calibrate_container, - update_calibrate_options=unmapped(dict(minuv=minuv)), + container=bandpass_options.calibrate_container, + update_calibrate_options=unmapped(dict(minuv=bandpass_options.minuv)), ) - for i in range(flag_calibrate_rounds): + for i in range(bandpass_options.flag_calibrate_rounds): # Apply and then recalibrate apply_cmds = task_bandpass_create_apply_solutions_cmd.map( ms=calibrate_cmds, calibrate_cmd=calibrate_cmds, output_column="CORRECTED_DATA", - container=calibrate_container, + container=bandpass_options.calibrate_container, ) flag_bandpass_mss = task_flag_ms_aoflagger.map( - ms=apply_cmds, container=flagger_container, rounds=1 + ms=apply_cmds, container=bandpass_options.flagger_container, rounds=1 ) calibrate_cmds = task_create_calibrate_cmd.map( ms=flag_bandpass_mss, calibrate_model=model_path, - container=calibrate_container, + container=bandpass_options.calibrate_container, calibrate_data_column="DATA", - update_calibrate_options=unmapped(dict(minuv=minuv)), + update_calibrate_options=unmapped(dict(minuv=bandpass_options.minuv)), ) flag_calibrate_cmds = task_flag_solutions.map( calibrate_cmd=calibrate_cmds, - smooth_window_size=smooth_window_size, - smooth_polynomial_order=smooth_polynomial_order, + smooth_window_size=bandpass_options.smooth_window_size, + smooth_polynomial_order=bandpass_options.smooth_polynomial_order, + mean_ant_tolerance=bandpass_options.mean_ant_tolerance, ) return flag_calibrate_cmds @@ -269,15 +262,10 @@ def calibrate_bandpass_flow( run_bandpass_stage( bandpass_mss=bandpass_mss, output_split_bandpass_path=output_split_bandpass_path, - calibrate_container=bandpass_options.calibrate_container, - flagger_container=bandpass_options.flagger_container, + bandpass_options=bandpass_options, model_path=model_path, source_name_prefix=source_name_prefix, skip_rotation=True, - smooth_window_size=bandpass_options.smooth_window_size, - smooth_polynomial_order=bandpass_options.smooth_polynomial_order, - flag_calibrate_rounds=bandpass_options.flag_calibrate_rounds, - minuv=bandpass_options.minuv, ) return output_split_bandpass_path @@ -399,6 +387,18 @@ def get_parser() -> ArgumentParser: default=None, help="The minimum baseline length, in meters, for data to be included in bandpass calibration stage", ) + parser.add_argument( + "--preflagger-ant-mean-tolerance", + type=float, + default=0.2, + help="Tolerance of the mean x/y antenna gain ratio test before antenna is flagged", + ) + parser.add_argument( + "--preflagger-mesh-ant-flags", + default=False, + action="store_true", + help="Share channel flags from bandpass solutions between all antennas", + ) return parser @@ -420,6 +420,8 @@ def cli() -> None: smooth_polynomial_order=args.smooth_polynomial_order, flag_calibrate_rounds=args.flag_calibrate_rounds, minuv=args.minuv, + preflagger_ant_mean_tolerance=args.preflagger_ant_mean_tolerance, + preflagger_mesh_ant_flags=args.preflagger_mesh_ant_flags, ) setup_run_bandpass_flow( diff --git a/tests/test_aocalibrate.py b/tests/test_aocalibrate.py index 7019700e..21e6f3cd 100644 --- a/tests/test_aocalibrate.py +++ b/tests/test_aocalibrate.py @@ -142,6 +142,27 @@ def test_sols_same_with_plots(ao_sols_known_bad): assert np.allclose(a_loaded.bandpass, b_loaded.bandpass, equal_nan=True) +def test_flagged_aosols_mesh(ao_sols_known_bad): + flagged_sols = flag_aosolutions( + solutions_path=ao_sols_known_bad, + plot_solutions_throughout=True, + smooth_solutions=True, + ) + assert isinstance(flagged_sols, FlaggedAOSolution) + assert len(flagged_sols.plots) == 9 + assert isinstance(flagged_sols.path, Path) + + flagged_sols = flag_aosolutions( + solutions_path=ao_sols_known_bad, + mesh_ant_flags=True, + plot_solutions_throughout=False, + smooth_solutions=False, + ) + assert isinstance(flagged_sols, FlaggedAOSolution) + assert len(flagged_sols.plots) == 0 + assert isinstance(flagged_sols.path, Path) + + def test_flagged_aosols(ao_sols_known_bad): flagged_sols = flag_aosolutions( solutions_path=ao_sols_known_bad, From 48c09bc493708dbebfe2ae73332d03c77d576963 Mon Sep 17 00:00:00 2001 From: tgalvin Date: Thu, 29 Feb 2024 16:01:10 +0800 Subject: [PATCH 02/21] corrected preflagger attribute name --- flint/prefect/flows/bandpass_pipeline.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/flint/prefect/flows/bandpass_pipeline.py b/flint/prefect/flows/bandpass_pipeline.py index 30fb3d17..779bc848 100644 --- a/flint/prefect/flows/bandpass_pipeline.py +++ b/flint/prefect/flows/bandpass_pipeline.py @@ -197,7 +197,8 @@ def run_bandpass_stage( calibrate_cmd=calibrate_cmds, smooth_window_size=bandpass_options.smooth_window_size, smooth_polynomial_order=bandpass_options.smooth_polynomial_order, - mean_ant_tolerance=bandpass_options.mean_ant_tolerance, + mean_ant_tolerance=bandpass_options.preflagger_mesh_ant_flagsmean_ant_tolerance, + mesh_ant_flags=bandpass_options.preflagger_mesh_ant_flags, ) return flag_calibrate_cmds From 27cc6e9c6a3b342858cce19e26492b3ba539d67d Mon Sep 17 00:00:00 2001 From: tgalvin Date: Thu, 29 Feb 2024 16:04:40 +0800 Subject: [PATCH 03/21] actually fixed attribute error --- flint/prefect/flows/bandpass_pipeline.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/flint/prefect/flows/bandpass_pipeline.py b/flint/prefect/flows/bandpass_pipeline.py index 779bc848..d7549282 100644 --- a/flint/prefect/flows/bandpass_pipeline.py +++ b/flint/prefect/flows/bandpass_pipeline.py @@ -197,7 +197,7 @@ def run_bandpass_stage( calibrate_cmd=calibrate_cmds, smooth_window_size=bandpass_options.smooth_window_size, smooth_polynomial_order=bandpass_options.smooth_polynomial_order, - mean_ant_tolerance=bandpass_options.preflagger_mesh_ant_flagsmean_ant_tolerance, + mean_ant_tolerance=bandpass_options.preflagger_mean_ant_tolerance, mesh_ant_flags=bandpass_options.preflagger_mesh_ant_flags, ) From 9f4a55c66805d291fdb06dad9c8b405b2b9dcdf0 Mon Sep 17 00:00:00 2001 From: tgalvin Date: Thu, 29 Feb 2024 16:13:41 +0800 Subject: [PATCH 04/21] no, now it has been fixed --- flint/prefect/flows/bandpass_pipeline.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/flint/prefect/flows/bandpass_pipeline.py b/flint/prefect/flows/bandpass_pipeline.py index d7549282..01de8af9 100644 --- a/flint/prefect/flows/bandpass_pipeline.py +++ b/flint/prefect/flows/bandpass_pipeline.py @@ -197,7 +197,7 @@ def run_bandpass_stage( calibrate_cmd=calibrate_cmds, smooth_window_size=bandpass_options.smooth_window_size, smooth_polynomial_order=bandpass_options.smooth_polynomial_order, - mean_ant_tolerance=bandpass_options.preflagger_mean_ant_tolerance, + mean_ant_tolerance=bandpass_options.preflagger_ant_mean_tolerance, mesh_ant_flags=bandpass_options.preflagger_mesh_ant_flags, ) From ff4797a517ab9dbb1f0fa183aa914d112b1f8dee Mon Sep 17 00:00:00 2001 From: tgalvin Date: Thu, 29 Feb 2024 16:22:36 +0800 Subject: [PATCH 05/21] typo in name --- flint/calibrate/aocalibrate.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/flint/calibrate/aocalibrate.py b/flint/calibrate/aocalibrate.py index 1ac63489..776c844b 100644 --- a/flint/calibrate/aocalibrate.py +++ b/flint/calibrate/aocalibrate.py @@ -864,7 +864,7 @@ def flag_aosolutions( plot_solutions_throughout: bool = True, smooth_window_size: int = 16, smooth_polynomial_order: int = 4, - mean_ant_tolderance: float = 0.1, + mean_ant_tolerance: float = 0.1, mesh_ant_flags: bool = False, ) -> FlaggedAOSolution: """Will open a previously solved ao-calibrate solutions file and flag additional channels and antennae. @@ -1002,7 +1002,7 @@ def flag_aosolutions( if flag_mean_xxyy_amplitude_ratio( xx_complex_gains=ant_gains[:, 0], yy_complex_gains=ant_gains[:, 3], - tolerance=mean_ant_tolderance, + tolerance=mean_ant_tolerance, ): logger.info(f"{ant=} failed mean amplitude gain test. Flagging {ant=}.") bandpass[time, ant, :, :] = np.nan From 7955580413b61808f6efe7709c7b3ff10121cacb Mon Sep 17 00:00:00 2001 From: tgalvin Date: Thu, 29 Feb 2024 16:27:05 +0800 Subject: [PATCH 06/21] passing through kwargs --- flint/prefect/flows/bandpass_pipeline.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/flint/prefect/flows/bandpass_pipeline.py b/flint/prefect/flows/bandpass_pipeline.py index 01de8af9..6fdbb1a5 100644 --- a/flint/prefect/flows/bandpass_pipeline.py +++ b/flint/prefect/flows/bandpass_pipeline.py @@ -81,6 +81,7 @@ def task_flag_solutions( calibrate_cmd: CalibrateCommand, smooth_window_size: int = 16, smooth_polynomial_order: int = 4, + **kwargs, ) -> CalibrateCommand: """Flag calibration solutions @@ -113,6 +114,7 @@ def task_flag_solutions( smooth_solutions=True, smooth_window_size=smooth_window_size, smooth_polynomial_order=smooth_polynomial_order, + **kwargs, ) for image_path in flagged_solutions.plots: From cabb95bfa2d8de79f4d4b3d80b61bd997d740d03 Mon Sep 17 00:00:00 2001 From: tgalvin Date: Thu, 29 Feb 2024 16:37:16 +0800 Subject: [PATCH 07/21] added some log messages --- flint/calibrate/aocalibrate.py | 8 +++++++- 1 file changed, 7 insertions(+), 1 deletion(-) diff --git a/flint/calibrate/aocalibrate.py b/flint/calibrate/aocalibrate.py index 776c844b..19cda79d 100644 --- a/flint/calibrate/aocalibrate.py +++ b/flint/calibrate/aocalibrate.py @@ -928,7 +928,13 @@ def flag_aosolutions( logger.info("Combining antenna flags") for time in range(solutions.nsol): time_mask = construct_mesh_ant_flags(mask=~np.isfinite(bandpass[time])) - bandpass[time, time_mask] = np.nan + logger.info( + f"Flags before applying mesh ant mask: {np.sum(~np.isfinite(bandpass[time]))}" + ) + bandpass[time][time_mask] = np.nan + logger.info( + f"Flags after applying mesh ant mask: {np.sum(~np.isfinite(bandpass[time]))}" + ) for time in range(solutions.nsol): for pol in (0, 3): From fb4f1f21177481818b4b11dcf23a3bcd08e6a935 Mon Sep 17 00:00:00 2001 From: tgalvin Date: Thu, 29 Feb 2024 16:38:00 +0800 Subject: [PATCH 08/21] default added --- flint/calibrate/aocalibrate.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/flint/calibrate/aocalibrate.py b/flint/calibrate/aocalibrate.py index 19cda79d..673b88ef 100644 --- a/flint/calibrate/aocalibrate.py +++ b/flint/calibrate/aocalibrate.py @@ -864,7 +864,7 @@ def flag_aosolutions( plot_solutions_throughout: bool = True, smooth_window_size: int = 16, smooth_polynomial_order: int = 4, - mean_ant_tolerance: float = 0.1, + mean_ant_tolerance: float = 0.2, mesh_ant_flags: bool = False, ) -> FlaggedAOSolution: """Will open a previously solved ao-calibrate solutions file and flag additional channels and antennae. @@ -887,7 +887,7 @@ def flag_aosolutions( plot_solutions_throughout (bool, Optional): If True, the solutions will be plotted at different stages of processing. Defaults to True. smooth_window_size (int, optional): The size of the window function of the savgol filter. Passed directly to savgol. Defaults to 16. smooth_polynomial_order (int, optional): The order of the polynomial of the savgol filter. Passed directly to savgol. Defaults to 4. - mean_ant_tolerance (float, optional): Tolerance of the mean x/y antenna gain ratio test before the antenna is flagged. Defaults to 0.1. + mean_ant_tolerance (float, optional): Tolerance of the mean x/y antenna gain ratio test before the antenna is flagged. Defaults to 0.2. mesh_ant_flags (bool, optional): If True, a channel is flagged across all antenna if it is flagged for any antenna. Defaults to False. Returns: From 3d2d60539facffc0b9dd75b1f96ba5a7ea03401e Mon Sep 17 00:00:00 2001 From: tgalvin Date: Thu, 29 Feb 2024 17:24:47 +0800 Subject: [PATCH 09/21] added tests --- flint/calibrate/aocalibrate.py | 2 +- tests/test_bptools.py | 35 ++++++++++++++++++++++++++++++++++ 2 files changed, 36 insertions(+), 1 deletion(-) create mode 100644 tests/test_bptools.py diff --git a/flint/calibrate/aocalibrate.py b/flint/calibrate/aocalibrate.py index 673b88ef..25271e3d 100644 --- a/flint/calibrate/aocalibrate.py +++ b/flint/calibrate/aocalibrate.py @@ -931,7 +931,7 @@ def flag_aosolutions( logger.info( f"Flags before applying mesh ant mask: {np.sum(~np.isfinite(bandpass[time]))}" ) - bandpass[time][time_mask] = np.nan + bandpass[time, time_mask] = np.nan logger.info( f"Flags after applying mesh ant mask: {np.sum(~np.isfinite(bandpass[time]))}" ) diff --git a/tests/test_bptools.py b/tests/test_bptools.py new file mode 100644 index 00000000..f48b4892 --- /dev/null +++ b/tests/test_bptools.py @@ -0,0 +1,35 @@ +"""Itemss around testing components of bptools +""" + +import pytest +import numpy as np + +from flint.bptools.preflagger import construct_mesh_ant_flags + + +def count_nan(data): + return np.sum(~np.isfinite(data)) + + +def test_construct_mesh_ant_flags(): + shape = (1, 36, 288, 4) + a = np.arange(np.prod(shape)).reshape(shape) * 1.0 + + assert count_nan(a) == 0 + + a[0, 0, 10:20, :] = np.nan + assert count_nan(a) == 40 + + mask = construct_mesh_ant_flags(mask=~np.isfinite(a[0])) + assert np.sum(mask) == 1440 + + a[0, mask] = np.nan + assert count_nan(a) == 1440 + + +def test_construct_mesh_ant_flags_assert(): + shape = (1, 36, 288, 4) + a = np.arange(np.prod(shape)).reshape(shape) * 1.0 + + with pytest.raises(AssertionError): + construct_mesh_ant_flags(mask=~np.isfinite(a)) From 4bbd5a3da57e1a3338b9528e7f687a1966b93285 Mon Sep 17 00:00:00 2001 From: tgalvin Date: Fri, 1 Mar 2024 10:03:15 +0800 Subject: [PATCH 10/21] added new flag and test --- flint/bptools/preflagger.py | 45 +++++++++++++++++++++++++++++++++ tests/test_bptools.py | 50 ++++++++++++++++++++++++++++++++++++- 2 files changed, 94 insertions(+), 1 deletion(-) diff --git a/flint/bptools/preflagger.py b/flint/bptools/preflagger.py index 39336c89..71234507 100644 --- a/flint/bptools/preflagger.py +++ b/flint/bptools/preflagger.py @@ -519,6 +519,8 @@ def construct_mesh_ant_flags(mask: np.ndarray) -> np.ndarray: empty_ants: List[int] = [] + # TODO: This can be replaced with numpy broadcasting + for ant in range(nant): ant_mask = mask[ant] if np.all(ant_mask): @@ -536,3 +538,46 @@ def construct_mesh_ant_flags(mask: np.ndarray) -> np.ndarray: result_mask[ant] = accumulate_mask return result_mask + + +def construct_jones_over_max_amp_flags( + complex_gains: np.ndarray, max_amplitude: float +) -> np.ndarray: + """Construct and return a mask that would flag an entire Jones + should there be an element whose amplitude is above a flagging + threshold + + Args: + complex_gains (np.ndarray): Complex gains that will have a mask constructed + max_amplitude (float): The flagging threshold, any Jones with a member above this will be flagged + + Returns: + np.ndarray: Boolean array of equal shape to `complex_gains`, with `True` indicating a flag + """ + + assert ( + complex_gains.shape[-1] == 4 + ), f"Expected last dimension to be length 4, received {complex_gains.shape=}" + + logger.info(f"Creating mask for Jones with amplitudes of {max_amplitude=}") + complex_gains = complex_gains.copy() + + original_shape = complex_gains.shape + + # Calculate tehe amplitudes of each of the complex numbers + # and construct the initial mask + amplitudes = np.abs(complex_gains) + mask = amplitudes > max_amplitude + + # Compress all but the last dimension into a single + # rank so we can easily broadcast over + mask = mask.reshape((-1, 4)) + + # Now broadcast like a pirate + flag_jones = np.any(mask, axis=1) + mask[flag_jones, :] = True + + # Convert back to original shape + mask = mask.reshape(original_shape) + + return mask diff --git a/tests/test_bptools.py b/tests/test_bptools.py index f48b4892..b1e60f12 100644 --- a/tests/test_bptools.py +++ b/tests/test_bptools.py @@ -4,7 +4,10 @@ import pytest import numpy as np -from flint.bptools.preflagger import construct_mesh_ant_flags +from flint.bptools.preflagger import ( + construct_mesh_ant_flags, + construct_jones_over_max_amp_flags, +) def count_nan(data): @@ -33,3 +36,48 @@ def test_construct_mesh_ant_flags_assert(): with pytest.raises(AssertionError): construct_mesh_ant_flags(mask=~np.isfinite(a)) + + +def test_construct_jobnes_over_max_amp_flags(): + shape = (1, 36, 288, 4) + a = np.arange(np.prod(shape)).reshape(shape) * 1.0 + + a[0, 10, 10, 0] = 100000000 + max_amplitude = 100000 + assert np.sum(a > max_amplitude) == 1 + + mask = construct_jones_over_max_amp_flags( + complex_gains=a, max_amplitude=max_amplitude + ) + assert np.sum(mask) == 4 + + a[mask] = np.nan + assert count_nan(a) == 4 + + +def test_construct_jobnes_over_max_amp_flags2(): + shape = (1, 36, 288, 4) + a = np.arange(np.prod(shape)).reshape(shape) * 1.0 + + a[0, 10, 10, 0] = 100000000 + a[0, 10, 11, 1] = 100000000 + a[0, 11, 12, 2] = 100000000 + a[0, 12, 13, 3] = 100000000 + max_amplitude = 100000 + assert np.sum(a > max_amplitude) == 4 + + mask = construct_jones_over_max_amp_flags( + complex_gains=a, max_amplitude=max_amplitude + ) + assert np.sum(mask) == 16 + + a[mask] = np.nan + assert count_nan(a) == 16 + + +def test_construct_jobnes_over_max_amp_flags3(): + shape = (1, 36, 288, 5) + a = np.arange(np.prod(shape)).reshape(shape) * 1.0 + + with pytest.raises(AssertionError): + construct_jones_over_max_amp_flags(complex_gains=a, max_amplitude=1000000) From 37c514ca687962819cc73d2c523d841b6153fbc8 Mon Sep 17 00:00:00 2001 From: tgalvin Date: Fri, 1 Mar 2024 10:18:54 +0800 Subject: [PATCH 11/21] added max flag mode to cli --- flint/calibrate/aocalibrate.py | 9 +++++++++ flint/options.py | 2 ++ flint/prefect/flows/bandpass_pipeline.py | 7 +++++++ 3 files changed, 18 insertions(+) diff --git a/flint/calibrate/aocalibrate.py b/flint/calibrate/aocalibrate.py index 25271e3d..08f61c51 100644 --- a/flint/calibrate/aocalibrate.py +++ b/flint/calibrate/aocalibrate.py @@ -27,6 +27,7 @@ flag_outlier_phase, flags_over_threshold, construct_mesh_ant_flags, + construct_jones_over_max_amp_flags, ) from flint.bptools.smoother import ( divide_bandpass_by_ref_ant_preserve_phase, @@ -866,6 +867,7 @@ def flag_aosolutions( smooth_polynomial_order: int = 4, mean_ant_tolerance: float = 0.2, mesh_ant_flags: bool = False, + max_gain_amplitude: Optional[float] = None, ) -> FlaggedAOSolution: """Will open a previously solved ao-calibrate solutions file and flag additional channels and antennae. @@ -889,6 +891,7 @@ def flag_aosolutions( smooth_polynomial_order (int, optional): The order of the polynomial of the savgol filter. Passed directly to savgol. Defaults to 4. mean_ant_tolerance (float, optional): Tolerance of the mean x/y antenna gain ratio test before the antenna is flagged. Defaults to 0.2. mesh_ant_flags (bool, optional): If True, a channel is flagged across all antenna if it is flagged for any antenna. Defaults to False. + max_gain_amplitude (Optional[float], optional): If not None, flag the Jones if an antenna has a amplitude gain above this value. Defaults to 10. Returns: FlaggedAOSolution: Path to the updated solutions file, intermediate solution files and plots along the way @@ -936,6 +939,12 @@ def flag_aosolutions( f"Flags after applying mesh ant mask: {np.sum(~np.isfinite(bandpass[time]))}" ) + if max_gain_amplitude: + mask = construct_jones_over_max_amp_flags( + complex_gains=bandpass, max_amplitude=max_gain_amplitude + ) + bandpass[mask] = np.nan + for time in range(solutions.nsol): for pol in (0, 3): logger.info(f"Processing {pols[pol]} polarisation") diff --git a/flint/options.py b/flint/options.py index 505e530b..6b4cdbab 100644 --- a/flint/options.py +++ b/flint/options.py @@ -37,6 +37,8 @@ class BandpassOptions(NamedTuple): """Tolerance that the mean x/y antenna gain ratio test before the antenna is flagged""" preflagger_mesh_ant_flags: bool = False """Share channel flags from bandpass solutions between all antenna""" + preflagger_jones_max_amplitudes: Optional[float] = None + """Flag Jones matrix if any amplitudes with a Jones are above this value""" class FieldOptions(NamedTuple): diff --git a/flint/prefect/flows/bandpass_pipeline.py b/flint/prefect/flows/bandpass_pipeline.py index 6fdbb1a5..83cdfce8 100644 --- a/flint/prefect/flows/bandpass_pipeline.py +++ b/flint/prefect/flows/bandpass_pipeline.py @@ -201,6 +201,7 @@ def run_bandpass_stage( smooth_polynomial_order=bandpass_options.smooth_polynomial_order, mean_ant_tolerance=bandpass_options.preflagger_ant_mean_tolerance, mesh_ant_flags=bandpass_options.preflagger_mesh_ant_flags, + max_fain_amplitude=bandpass_options.preflagger_jones_max_amplitude, ) return flag_calibrate_cmds @@ -402,6 +403,12 @@ def get_parser() -> ArgumentParser: action="store_true", help="Share channel flags from bandpass solutions between all antennas", ) + parser.add_argument( + "--preflagger-jones-max-amplitude", + default=None, + type=float, + help="Flag Jones matrix if any amplitudes with a Jones are above this value", + ) return parser From 454b74e771b1c64f6ed82bbe2e08144018d0ffc4 Mon Sep 17 00:00:00 2001 From: tgalvin Date: Fri, 1 Mar 2024 10:34:08 +0800 Subject: [PATCH 12/21] updated option name --- flint/options.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/flint/options.py b/flint/options.py index 6b4cdbab..0cef03e4 100644 --- a/flint/options.py +++ b/flint/options.py @@ -37,7 +37,7 @@ class BandpassOptions(NamedTuple): """Tolerance that the mean x/y antenna gain ratio test before the antenna is flagged""" preflagger_mesh_ant_flags: bool = False """Share channel flags from bandpass solutions between all antenna""" - preflagger_jones_max_amplitudes: Optional[float] = None + preflagger_jones_max_amplitude: Optional[float] = None """Flag Jones matrix if any amplitudes with a Jones are above this value""" From 7ddb6d5628c86ef80e265ad0fd5a942275434986 Mon Sep 17 00:00:00 2001 From: tgalvin Date: Fri, 1 Mar 2024 10:38:59 +0800 Subject: [PATCH 13/21] corrected typo --- flint/prefect/flows/bandpass_pipeline.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/flint/prefect/flows/bandpass_pipeline.py b/flint/prefect/flows/bandpass_pipeline.py index 83cdfce8..503ac57c 100644 --- a/flint/prefect/flows/bandpass_pipeline.py +++ b/flint/prefect/flows/bandpass_pipeline.py @@ -201,7 +201,7 @@ def run_bandpass_stage( smooth_polynomial_order=bandpass_options.smooth_polynomial_order, mean_ant_tolerance=bandpass_options.preflagger_ant_mean_tolerance, mesh_ant_flags=bandpass_options.preflagger_mesh_ant_flags, - max_fain_amplitude=bandpass_options.preflagger_jones_max_amplitude, + max_gain_amplitude=bandpass_options.preflagger_jones_max_amplitude, ) return flag_calibrate_cmds From db4311aefaa9de19678dad01e36bee883e816933 Mon Sep 17 00:00:00 2001 From: tgalvin Date: Fri, 1 Mar 2024 10:44:43 +0800 Subject: [PATCH 14/21] added logger message when flagging --- flint/bptools/preflagger.py | 5 +++++ flint/calibrate/aocalibrate.py | 1 + 2 files changed, 6 insertions(+) diff --git a/flint/bptools/preflagger.py b/flint/bptools/preflagger.py index 71234507..489b491f 100644 --- a/flint/bptools/preflagger.py +++ b/flint/bptools/preflagger.py @@ -580,4 +580,9 @@ def construct_jones_over_max_amp_flags( # Convert back to original shape mask = mask.reshape(original_shape) + no_flagged = np.sum(mask) + + if no_flagged > 0: + logger.warning(f"{no_flagged} items flagged with {max_amplitude=}") + return mask diff --git a/flint/calibrate/aocalibrate.py b/flint/calibrate/aocalibrate.py index 08f61c51..a7bf7b12 100644 --- a/flint/calibrate/aocalibrate.py +++ b/flint/calibrate/aocalibrate.py @@ -943,6 +943,7 @@ def flag_aosolutions( mask = construct_jones_over_max_amp_flags( complex_gains=bandpass, max_amplitude=max_gain_amplitude ) + no_flagged = np.sum(mask) bandpass[mask] = np.nan for time in range(solutions.nsol): From ee7c238db7ca3063d3d1b257b1bd8cf60fc5ba77 Mon Sep 17 00:00:00 2001 From: tgalvin Date: Fri, 1 Mar 2024 11:10:13 +0800 Subject: [PATCH 15/21] changed accumulate mask slicing --- flint/bptools/preflagger.py | 2 ++ flint/calibrate/aocalibrate.py | 16 ++++++---------- 2 files changed, 8 insertions(+), 10 deletions(-) diff --git a/flint/bptools/preflagger.py b/flint/bptools/preflagger.py index 489b491f..d963e153 100644 --- a/flint/bptools/preflagger.py +++ b/flint/bptools/preflagger.py @@ -529,6 +529,8 @@ def construct_mesh_ant_flags(mask: np.ndarray) -> np.ndarray: accumulate_mask = accumulate_mask | ant_mask + logger.info(f"Flags in accumulated mask: {np.sum(accumulate_mask)}") + result_mask = np.zeros_like(mask, dtype=bool) for ant in range(nant): diff --git a/flint/calibrate/aocalibrate.py b/flint/calibrate/aocalibrate.py index a7bf7b12..57882563 100644 --- a/flint/calibrate/aocalibrate.py +++ b/flint/calibrate/aocalibrate.py @@ -857,7 +857,7 @@ class FlaggedAOSolution(NamedTuple): def flag_aosolutions( solutions_path: Path, - ref_ant: Optional[int] = -1, + ref_ant: int = -1, flag_cut: float = 3, plot_dir: Optional[Path] = None, out_solutions_path: Optional[Path] = None, @@ -929,21 +929,17 @@ def flag_aosolutions( if mesh_ant_flags: logger.info("Combining antenna flags") + mask = np.zeros_like(bandpass, dtype=bool) + for time in range(solutions.nsol): - time_mask = construct_mesh_ant_flags(mask=~np.isfinite(bandpass[time])) - logger.info( - f"Flags before applying mesh ant mask: {np.sum(~np.isfinite(bandpass[time]))}" - ) - bandpass[time, time_mask] = np.nan - logger.info( - f"Flags after applying mesh ant mask: {np.sum(~np.isfinite(bandpass[time]))}" - ) + mask[time] = construct_mesh_ant_flags(mask=~np.isfinite(bandpass[time])) + + bandpass[mask] = np.nan if max_gain_amplitude: mask = construct_jones_over_max_amp_flags( complex_gains=bandpass, max_amplitude=max_gain_amplitude ) - no_flagged = np.sum(mask) bandpass[mask] = np.nan for time in range(solutions.nsol): From aaab11fe3493e42426c7adfe902b92eaa3d4cc49 Mon Sep 17 00:00:00 2001 From: tgalvin Date: Fri, 1 Mar 2024 11:19:43 +0800 Subject: [PATCH 16/21] added --smooth-solutions --- flint/options.py | 2 ++ flint/prefect/flows/bandpass_pipeline.py | 15 ++++++++------- 2 files changed, 10 insertions(+), 7 deletions(-) diff --git a/flint/options.py b/flint/options.py index 0cef03e4..ead2b187 100644 --- a/flint/options.py +++ b/flint/options.py @@ -25,6 +25,8 @@ class BandpassOptions(NamedTuple): """Path to the singularity calibrate container""" expected_ms: int = 36 """The expected number of measurement set files to find""" + smooth_solutions: bool = False + """Will activate the smoothing of the bandpass solutions""" smooth_window_size: int = 16 """The width of the smoothing window used to smooth the bandpass solutions""" smooth_polynomial_order: int = 4 diff --git a/flint/prefect/flows/bandpass_pipeline.py b/flint/prefect/flows/bandpass_pipeline.py index 503ac57c..4d1b56a4 100644 --- a/flint/prefect/flows/bandpass_pipeline.py +++ b/flint/prefect/flows/bandpass_pipeline.py @@ -79,16 +79,12 @@ def task_bandpass_create_apply_solutions_cmd( @task def task_flag_solutions( calibrate_cmd: CalibrateCommand, - smooth_window_size: int = 16, - smooth_polynomial_order: int = 4, **kwargs, ) -> CalibrateCommand: """Flag calibration solutions Args: calibrate_cmd (CalibrateCommand): Calibrate command that contains path to the solution file that will be flagged - smooth_window_size (int, optional): The size of the window function of the savgol filter. Passed directly to savgol. Defaults to 16. - smooth_polynomial_order (int, optional): The order of the polynomial of the savgol filter. Passed directly to savgol. Defaults to 4. Returns: CalibrateCommand: Calibrate command with update meta-data describing the new solutions file @@ -111,9 +107,6 @@ def task_flag_solutions( ref_ant=-1, flag_cut=3, plot_dir=plot_dir, - smooth_solutions=True, - smooth_window_size=smooth_window_size, - smooth_polynomial_order=smooth_polynomial_order, **kwargs, ) @@ -197,6 +190,7 @@ def run_bandpass_stage( ) flag_calibrate_cmds = task_flag_solutions.map( calibrate_cmd=calibrate_cmds, + smooth_solutions=bandpass_options.smooth_solutions, smooth_window_size=bandpass_options.smooth_window_size, smooth_polynomial_order=bandpass_options.smooth_polynomial_order, mean_ant_tolerance=bandpass_options.preflagger_ant_mean_tolerance, @@ -367,6 +361,12 @@ def get_parser() -> ArgumentParser: default="petrichor", help="Path to a cluster configuration file, or a known cluster name. ", ) + parser.add_argument( + "--smooth-solutions", + default=False, + action="store_true", + help="Smooth the bandpass solutions", + ) parser.add_argument( "--smooth-window-size", default=16, @@ -426,6 +426,7 @@ def cli() -> None: flagger_container=args.flagger_container, calibrate_container=args.calibrate_container, expected_ms=args.expected_ms, + smooth_solutions=args.smooth_solutions, smooth_window_size=args.smooth_window_size, smooth_polynomial_order=args.smooth_polynomial_order, flag_calibrate_rounds=args.flag_calibrate_rounds, From 03a25c69fa38e98b50a1de5ec9460658401ca42c Mon Sep 17 00:00:00 2001 From: tgalvin Date: Fri, 1 Mar 2024 14:37:01 +0800 Subject: [PATCH 17/21] added some notes and docstring --- flint/calibrate/aocalibrate.py | 25 ++++++++++++++++++++++--- 1 file changed, 22 insertions(+), 3 deletions(-) diff --git a/flint/calibrate/aocalibrate.py b/flint/calibrate/aocalibrate.py index 57882563..1d36d2aa 100644 --- a/flint/calibrate/aocalibrate.py +++ b/flint/calibrate/aocalibrate.py @@ -871,11 +871,28 @@ def flag_aosolutions( ) -> FlaggedAOSolution: """Will open a previously solved ao-calibrate solutions file and flag additional channels and antennae. - There are currently two main stages. The first will attempt to search for channels where the the phase of the + There are a number of distinct operations applied to the data, which are + presented in order they are applied. + + If `mesh_ant_flags` is `True`, channels flagged from on channel on a single + antenna will be applied to all (unless an antenna is completely flagged). + This happens before any other operation,. + + If `max_gain_amplitude` is not `None` than any Jones with an element + whose amplitude is above the set value will be flagged. + + Next, an attempt is made to search for channels where the the phase of the gain solution are outliers. The phase over frequency is first unwrapped (delay solved for) before the flagging statistics are computed. - The second stage will flag an entire antenna if more then 80 percent of the flags for a polarisation are flagged. + If an antenna is over 80% flagged then it is completely removed. + + A low order polynomial (typically order 5) is fit to the amplitudes of the + Gx and Gy, and if the residuals are sufficently high then the pol will + be flagged. + + If the mean ratio of the Gx and Gy amplitudes for an antenna are higher + then `mean_ant_tolerance` then the antenna will be flagged. Keywords that with the `smooth` prefix are passed to the `smooth_bandpass_complex_gains` function. @@ -890,13 +907,15 @@ def flag_aosolutions( smooth_window_size (int, optional): The size of the window function of the savgol filter. Passed directly to savgol. Defaults to 16. smooth_polynomial_order (int, optional): The order of the polynomial of the savgol filter. Passed directly to savgol. Defaults to 4. mean_ant_tolerance (float, optional): Tolerance of the mean x/y antenna gain ratio test before the antenna is flagged. Defaults to 0.2. - mesh_ant_flags (bool, optional): If True, a channel is flagged across all antenna if it is flagged for any antenna. Defaults to False. + mesh_ant_flags (bool, optional): If True, a channel is flagged across all antenna if it is flagged for any antenna. Performed before other flagging operations. Defaults to False. max_gain_amplitude (Optional[float], optional): If not None, flag the Jones if an antenna has a amplitude gain above this value. Defaults to 10. Returns: FlaggedAOSolution: Path to the updated solutions file, intermediate solution files and plots along the way """ # TODO: This should be broken down into separate stages. Way too large of a function. + # TODO: This pirate needs to cull some of this logic out, likely not needed + # and dead solutions = AOSolutions.load(path=solutions_path) title = solutions_path.name From c725041fe1d9687af509cfdb14531710d0764683 Mon Sep 17 00:00:00 2001 From: tgalvin Date: Fri, 1 Mar 2024 14:55:29 +0800 Subject: [PATCH 18/21] added helper create_directory function, slight clean --- flint/calibrate/aocalibrate.py | 22 +++++++++++----------- flint/utils.py | 27 +++++++++++++++++++++++++++ 2 files changed, 38 insertions(+), 11 deletions(-) diff --git a/flint/calibrate/aocalibrate.py b/flint/calibrate/aocalibrate.py index 1d36d2aa..a142e5a5 100644 --- a/flint/calibrate/aocalibrate.py +++ b/flint/calibrate/aocalibrate.py @@ -38,6 +38,7 @@ from flint.ms import MS, consistent_ms, get_beam_from_ms from flint.naming import get_aocalibrate_output_path from flint.sclient import run_singularity_command +from flint.utils import create_directory class CalibrateOptions(NamedTuple): @@ -888,7 +889,7 @@ def flag_aosolutions( If an antenna is over 80% flagged then it is completely removed. A low order polynomial (typically order 5) is fit to the amplitudes of the - Gx and Gy, and if the residuals are sufficently high then the pol will + Gx and Gy, and if the residuals are sufficently high then the antenna will be flagged. If the mean ratio of the Gx and Gy amplitudes for an antenna are higher @@ -922,12 +923,8 @@ def flag_aosolutions( pols = {0: "XX", 1: "XY", 2: "YX", 3: "YY"} - if plot_dir is not None and not plot_dir.exists(): - logger.info(f"Creating {str(plot_dir)}") - try: - plot_dir.mkdir(parents=True) - except Exception as e: - logger.error(f"Failed to create {str(plot_dir)} {e}.") + if plot_dir: + create_directory(directory=plot_dir) # Note that although the solutions variable (an instance of AOSolutions) is immutable, # which includes the reference to the numpy array, the _actual_ numpy array is! So, @@ -962,6 +959,9 @@ def flag_aosolutions( bandpass[mask] = np.nan for time in range(solutions.nsol): + ref_bandpass = divide_bandpass_by_ref_ant_preserve_phase( + complex_gains=bandpass[time], ref_ant=ref_ant + ) for pol in (0, 3): logger.info(f"Processing {pols[pol]} polarisation") ref_ant_gains = bandpass[time, ref_ant, :, pol] @@ -973,7 +973,7 @@ def flag_aosolutions( logger.info(f"Skipping reference antenna = ant{ref_ant:02}") continue - ant_gains = bandpass[time, ant, :, pol] / ref_ant_gains + ant_gains = ref_bandpass[ant] plot_title = f"{title} - ant{ant:02d} - {pols[pol]}" ouput_path = ( plot_dir / f"{title}.ant{ant:02d}.{pols[pol]}.png" @@ -995,7 +995,7 @@ def flag_aosolutions( bandpass[time, ant, phase_outlier_result.outlier_mask, pol] = np.nan except PhaseOutlierFitError: # This is raised if the fit failed to converge, or some other nasty. - bandpass[time, ant, :, pol] = np.nan + bandpass[time, ant, :, :] = np.nan # Flag all solutions for this (ant,pol) if more than 80% are flagged if flags_over_threshold( @@ -1006,7 +1006,7 @@ def flag_aosolutions( logger.info( f"Flagging all solutions across {pols[pol]} for ant{ant:02d}, too many flagged channels." ) - bandpass[time, ant, :, pol] = np.nan + bandpass[time, ant, :, :] = np.nan complex_gains = bandpass[time, ant, :, pol] if any(np.isfinite(complex_gains)) and flag_mean_residual_amplitude( @@ -1015,7 +1015,7 @@ def flag_aosolutions( logger.info( f"Flagging all solutions across {pols[pol]} for ant{ant:02d}, mean residual amplitudes high" ) - bandpass[time, ant, :, pol] = np.nan + bandpass[time, ant, :, :] = np.nan flagged = ~np.isfinite(bandpass[time, ant, :, pol]) logger.info( diff --git a/flint/utils.py b/flint/utils.py index 991417eb..990418c0 100644 --- a/flint/utils.py +++ b/flint/utils.py @@ -157,3 +157,30 @@ def remove_files_folders(*paths_to_remove: Path) -> List[Path]: files_removed.append(file) return files_removed + + +def create_directory(directory: Path) -> Path: + """Will attempt to safely create a directory. Should it + not exist it will be created. if this creates an exception, + which might happen in a multi-process environment, it is + reported and discarded. + + Args: + directory (Path): Path to directory to create + + Returns: + Path: The directory created + """ + + directory = Path(directory) + + if directory.exists(): + return directory + + logger.info(f"Creating {str(directory)}") + try: + directory.mkdir(parents=True) + except Exception as e: + logger.error(f"Failed to create {str(directory)} {e}.") + + return directory From a170d0bdfd524d592a2767ad5471a81089763f34 Mon Sep 17 00:00:00 2001 From: Tim Galvin Date: Fri, 1 Mar 2024 18:34:34 +1100 Subject: [PATCH 19/21] fixed overlooked phase outlier --- flint/bptools/preflagger.py | 3 +++ flint/calibrate/aocalibrate.py | 17 +++++++++-------- 2 files changed, 12 insertions(+), 8 deletions(-) diff --git a/flint/bptools/preflagger.py b/flint/bptools/preflagger.py index d963e153..73c267eb 100644 --- a/flint/bptools/preflagger.py +++ b/flint/bptools/preflagger.py @@ -405,6 +405,9 @@ def flag_mean_residual_amplitude( bool: Whether the data should be considered bad. True if it is bad, False if otherwise. """ + if not np.any(np.isfinite(complex_gains)): + return True + amplitudes = np.abs(complex_gains) idxs = np.arange(amplitudes.shape[0]) mask = np.isfinite(amplitudes) diff --git a/flint/calibrate/aocalibrate.py b/flint/calibrate/aocalibrate.py index a142e5a5..5f2e11f9 100644 --- a/flint/calibrate/aocalibrate.py +++ b/flint/calibrate/aocalibrate.py @@ -965,15 +965,13 @@ def flag_aosolutions( for pol in (0, 3): logger.info(f"Processing {pols[pol]} polarisation") ref_ant_gains = bandpass[time, ref_ant, :, pol] - if np.sum(np.isfinite(ref_ant_gains)) == 0: - raise ValueError(f"The ref_ant={ref_ant} is completely bad. ") - + for ant in range(solutions.nant): if ant == ref_ant: logger.info(f"Skipping reference antenna = ant{ref_ant:02}") continue - ant_gains = ref_bandpass[ant] + ant_gains = ref_bandpass[ant, :, pol] plot_title = f"{title} - ant{ant:02d} - {pols[pol]}" ouput_path = ( plot_dir / f"{title}.ant{ant:02d}.{pols[pol]}.png" @@ -992,11 +990,14 @@ def flag_aosolutions( plot_title=plot_title, plot_path=ouput_path, ) - bandpass[time, ant, phase_outlier_result.outlier_mask, pol] = np.nan + bandpass[time, ant, phase_outlier_result.outlier_mask, :] = np.nan except PhaseOutlierFitError: # This is raised if the fit failed to converge, or some other nasty. bandpass[time, ant, :, :] = np.nan + for time in range(solutions.nsol): + for pol in (0, 3): + for ant in range(solutions.nant): # Flag all solutions for this (ant,pol) if more than 80% are flagged if flags_over_threshold( flags=~np.isfinite(bandpass[time, ant, :, pol]), @@ -1004,16 +1005,16 @@ def flag_aosolutions( ant_idx=ant, ): logger.info( - f"Flagging all solutions across {pols[pol]} for ant{ant:02d}, too many flagged channels." + f"Flagging all solutions across ant{ant:02d}, too many flagged channels." ) bandpass[time, ant, :, :] = np.nan complex_gains = bandpass[time, ant, :, pol] - if any(np.isfinite(complex_gains)) and flag_mean_residual_amplitude( + if flag_mean_residual_amplitude( complex_gains=complex_gains ): logger.info( - f"Flagging all solutions across {pols[pol]} for ant{ant:02d}, mean residual amplitudes high" + f"Flagging all solutions for ant{ant:02d}, mean residual amplitudes high" ) bandpass[time, ant, :, :] = np.nan From 794b48b415e874153639ba78bcef5844afaa963b Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Fri, 1 Mar 2024 08:16:10 +0000 Subject: [PATCH 20/21] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- flint/bptools/preflagger.py | 2 +- flint/calibrate/aocalibrate.py | 11 ++++------- tests/test_bptools.py | 4 ++-- 3 files changed, 7 insertions(+), 10 deletions(-) diff --git a/flint/bptools/preflagger.py b/flint/bptools/preflagger.py index 73c267eb..8d1ddb5c 100644 --- a/flint/bptools/preflagger.py +++ b/flint/bptools/preflagger.py @@ -7,7 +7,7 @@ """ from pathlib import Path -from typing import NamedTuple, Optional, Tuple, List +from typing import List, NamedTuple, Optional, Tuple import matplotlib.pyplot as plt import numpy as np diff --git a/flint/calibrate/aocalibrate.py b/flint/calibrate/aocalibrate.py index 5f2e11f9..f63dfefc 100644 --- a/flint/calibrate/aocalibrate.py +++ b/flint/calibrate/aocalibrate.py @@ -22,12 +22,12 @@ import numpy as np from flint.bptools.preflagger import ( + construct_jones_over_max_amp_flags, + construct_mesh_ant_flags, flag_mean_residual_amplitude, flag_mean_xxyy_amplitude_ratio, flag_outlier_phase, flags_over_threshold, - construct_mesh_ant_flags, - construct_jones_over_max_amp_flags, ) from flint.bptools.smoother import ( divide_bandpass_by_ref_ant_preserve_phase, @@ -965,7 +965,7 @@ def flag_aosolutions( for pol in (0, 3): logger.info(f"Processing {pols[pol]} polarisation") ref_ant_gains = bandpass[time, ref_ant, :, pol] - + for ant in range(solutions.nant): if ant == ref_ant: logger.info(f"Skipping reference antenna = ant{ref_ant:02}") @@ -1010,9 +1010,7 @@ def flag_aosolutions( bandpass[time, ant, :, :] = np.nan complex_gains = bandpass[time, ant, :, pol] - if flag_mean_residual_amplitude( - complex_gains=complex_gains - ): + if flag_mean_residual_amplitude(complex_gains=complex_gains): logger.info( f"Flagging all solutions for ant{ant:02d}, mean residual amplitudes high" ) @@ -1029,7 +1027,6 @@ def flag_aosolutions( ) # This loop will flag based on stats across different polarisations for ant in range(solutions.nant): - ant_gains = bandpass_phased_referenced[ant] if flag_mean_xxyy_amplitude_ratio( xx_complex_gains=ant_gains[:, 0], diff --git a/tests/test_bptools.py b/tests/test_bptools.py index b1e60f12..c8e774e6 100644 --- a/tests/test_bptools.py +++ b/tests/test_bptools.py @@ -1,12 +1,12 @@ """Itemss around testing components of bptools """ -import pytest import numpy as np +import pytest from flint.bptools.preflagger import ( - construct_mesh_ant_flags, construct_jones_over_max_amp_flags, + construct_mesh_ant_flags, ) From 111502b4a041c977fd1fb31936cef4adb568cc36 Mon Sep 17 00:00:00 2001 From: Tim Galvin Date: Fri, 1 Mar 2024 19:27:57 +1100 Subject: [PATCH 21/21] removed unused variable --- flint/calibrate/aocalibrate.py | 1 - 1 file changed, 1 deletion(-) diff --git a/flint/calibrate/aocalibrate.py b/flint/calibrate/aocalibrate.py index 5f2e11f9..83584e04 100644 --- a/flint/calibrate/aocalibrate.py +++ b/flint/calibrate/aocalibrate.py @@ -964,7 +964,6 @@ def flag_aosolutions( ) for pol in (0, 3): logger.info(f"Processing {pols[pol]} polarisation") - ref_ant_gains = bandpass[time, ref_ant, :, pol] for ant in range(solutions.nant): if ant == ref_ant: