From ce1db3e9808cc97a3185a29f22d070c7b2bb601f Mon Sep 17 00:00:00 2001 From: "Thomson, Alec (CASS, Kensington)" Date: Mon, 13 May 2024 11:06:49 +1000 Subject: [PATCH] Improve boundary checker --- arrakis/cutout.py | 30 +++++++++++++++++++++++------- 1 file changed, 23 insertions(+), 7 deletions(-) diff --git a/arrakis/cutout.py b/arrakis/cutout.py index 712b8d49..5f9ec62c 100644 --- a/arrakis/cutout.py +++ b/arrakis/cutout.py @@ -256,9 +256,23 @@ def get_args( coords = SkyCoord(ras, decs) try: + # If in South - Northern boundary will have greatest RA range + # And vice versa northern_boundary = SkyCoord(ras, np.ones_like(ras.value) * decs.max()) - for i, coord_1 in enumerate(northern_boundary): - for j, coord_2 in enumerate(northern_boundary): + southern_boundary = SkyCoord(ras, np.ones_like(ras.value) * decs.min()) + + if abs(decs.max()) > abs(decs.min()): + # North - use Southern + longest_boundary = southern_boundary + else: + # South - use Northern + longest_boundary = northern_boundary + + # Find the coords with greatest separation + # That gives us our image size + # PA tells us East from West wrap + for i, coord_1 in enumerate(longest_boundary): + for j, coord_2 in enumerate(longest_boundary): sep = coord_1.separation(coord_2) pa = coord_1.position_angle(coord_2) if i == 0 and j == 0: @@ -272,12 +286,14 @@ def get_args( largest_sep = sep largest_pa = pa - if largest_pa > 90 * u.deg: - left_point = northern_boundary[largest_i] - right_point = northern_boundary[largest_j] + if largest_pa > 180 * u.deg: + # PA points West / left + left_point = longest_boundary[largest_i] + right_point = longest_boundary[largest_j] else: - left_point = northern_boundary[largest_j] - right_point = northern_boundary[largest_i] + # PA points East / right + left_point = longest_boundary[largest_j] + right_point = longest_boundary[largest_i] # Use SkyCoords to account for poles and meridian top_right = SkyCoord(right_point.ra, decs.max())