From af405a62c28a4af6792dc8cb01710024ad8e5741 Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Wed, 18 Jun 2025 20:40:28 +0100 Subject: [PATCH 01/19] resized_array_2d_from converted to numpy --- autoarray/operators/transformer.py | 8 +- autoarray/structures/arrays/array_2d_util.py | 73 +++++++------------ .../structures/arrays/test_array_2d_util.py | 1 - 3 files changed, 29 insertions(+), 53 deletions(-) diff --git a/autoarray/operators/transformer.py b/autoarray/operators/transformer.py index c3d94f686..675123e78 100644 --- a/autoarray/operators/transformer.py +++ b/autoarray/operators/transformer.py @@ -422,13 +422,13 @@ def transform_mapping_matrix(self, mapping_matrix: np.ndarray) -> np.ndarray: Parameters ---------- - mapping_matrix : np.ndarray + mapping_matrix A 2D array where each column corresponds to a source-plane pixel intensity distribution flattened into image space. Returns - ------- - np.ndarray - A complex-valued 2D array where each column contains the visibilities corresponding to the respective column in the input mapping matrix. + ------- + A complex-valued 2D array where each column contains the visibilities corresponding to the respective column + in the input mapping matrix. Notes ----- diff --git a/autoarray/structures/arrays/array_2d_util.py b/autoarray/structures/arrays/array_2d_util.py index 4efa6583b..c006e3d9f 100644 --- a/autoarray/structures/arrays/array_2d_util.py +++ b/autoarray/structures/arrays/array_2d_util.py @@ -272,11 +272,10 @@ def extracted_array_2d_from( return resized_array -@numba_util.jit() def resized_array_2d_from( array_2d: np.ndarray, resized_shape: Tuple[int, int], - origin: Tuple[int, int] = (-1, -1), + origin: Tuple[int, int] = None, pad_value: int = 0.0, ) -> np.ndarray: """ @@ -312,56 +311,34 @@ def resized_array_2d_from( resize_array = resize_array_2d(array_2d=array_2d, new_shape=(2,2), origin=(2, 2)) """ - y_is_even = int(array_2d.shape[0]) % 2 == 0 - x_is_even = int(array_2d.shape[1]) % 2 == 0 + if origin is None: + y_centre = array_2d.shape[0] // 2 + x_centre = array_2d.shape[1] // 2 + origin = (y_centre, x_centre) - if origin == (-1, -1): - if y_is_even: - y_centre = int(array_2d.shape[0] / 2) - elif not y_is_even: - y_centre = int(array_2d.shape[0] / 2) + # Define window edges so that length == resized_shape dimension exactly + y_min = origin[0] - resized_shape[0] // 2 + y_max = y_min + resized_shape[0] - if x_is_even: - x_centre = int(array_2d.shape[1] / 2) - elif not x_is_even: - x_centre = int(array_2d.shape[1] / 2) + x_min = origin[1] - resized_shape[1] // 2 + x_max = x_min + resized_shape[1] - origin = (y_centre, x_centre) + resized_array = np.full(resized_shape, pad_value, dtype=array_2d.dtype) + + # Calculate source indices clipped to array bounds + src_y_start = max(y_min, 0) + src_y_end = min(y_max, array_2d.shape[0]) + src_x_start = max(x_min, 0) + src_x_end = min(x_max, array_2d.shape[1]) + + # Calculate destination indices corresponding to source indices + dst_y_start = max(0, -y_min) + dst_y_end = dst_y_start + (src_y_end - src_y_start) + dst_x_start = max(0, -x_min) + dst_x_end = dst_x_start + (src_x_end - src_x_start) - resized_array = np.zeros(shape=resized_shape) - - if y_is_even: - y_min = origin[0] - int(resized_shape[0] / 2) - y_max = origin[0] + int((resized_shape[0] / 2)) + 1 - elif not y_is_even: - y_min = origin[0] - int(resized_shape[0] / 2) - y_max = origin[0] + int((resized_shape[0] / 2)) + 1 - - if x_is_even: - x_min = origin[1] - int(resized_shape[1] / 2) - x_max = origin[1] + int((resized_shape[1] / 2)) + 1 - elif not x_is_even: - x_min = origin[1] - int(resized_shape[1] / 2) - x_max = origin[1] + int((resized_shape[1] / 2)) + 1 - - for y_resized, y in enumerate(range(y_min, y_max)): - for x_resized, x in enumerate(range(x_min, x_max)): - if y >= 0 and y < array_2d.shape[0] and x >= 0 and x < array_2d.shape[1]: - if ( - y_resized >= 0 - and y_resized < resized_shape[0] - and x_resized >= 0 - and x_resized < resized_shape[1] - ): - resized_array[y_resized, x_resized] = array_2d[y, x] - else: - if ( - y_resized >= 0 - and y_resized < resized_shape[0] - and x_resized >= 0 - and x_resized < resized_shape[1] - ): - resized_array[y_resized, x_resized] = pad_value + # Copy overlapping region from source to destination + resized_array[dst_y_start:dst_y_end, dst_x_start:dst_x_end] = array_2d[src_y_start:src_y_end, src_x_start:src_x_end] return resized_array diff --git a/test_autoarray/structures/arrays/test_array_2d_util.py b/test_autoarray/structures/arrays/test_array_2d_util.py index 9468dd510..0fe0d0a01 100644 --- a/test_autoarray/structures/arrays/test_array_2d_util.py +++ b/test_autoarray/structures/arrays/test_array_2d_util.py @@ -247,7 +247,6 @@ def test__resized_array_2d_from__padding(): ) ).all() - def test__resized_array_2d_from__padding_with_new_origin(): array = np.ones((3, 3)) array[2, 2] = 2.0 From f095be5ea389dd11e60270f2fbbe76168b0bd7b3 Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Wed, 18 Jun 2025 20:41:15 +0100 Subject: [PATCH 02/19] remove replace_noise_map_2d_values_where_image_2d_values_are_negative --- autoarray/structures/arrays/array_2d_util.py | 45 ---------- .../structures/arrays/test_array_2d_util.py | 89 ------------------- 2 files changed, 134 deletions(-) diff --git a/autoarray/structures/arrays/array_2d_util.py b/autoarray/structures/arrays/array_2d_util.py index c006e3d9f..a9be13fe6 100644 --- a/autoarray/structures/arrays/array_2d_util.py +++ b/autoarray/structures/arrays/array_2d_util.py @@ -343,51 +343,6 @@ def resized_array_2d_from( return resized_array -@numba_util.jit() -def replace_noise_map_2d_values_where_image_2d_values_are_negative( - image_2d: np.ndarray, noise_map_2d: np.ndarray, target_signal_to_noise: float = 2.0 -) -> np.ndarray: - """ - If the values of a 2D image array are negative, this function replaces the corresponding 2D noise-map array - values to meet a specified target to noise value. - - This routine is necessary because of anomolous values in images which come from our HST ACS data_type-reduction - pipeline, where image-pixels with negative values (e.g. due to the background sky subtraction) have extremely - small noise values, which inflate their signal-to-noise values and chi-squared contributions in the modeling. - - Parameters - ---------- - image_2d - The 2D image array used to locate the pixel indexes in the noise-map which are replaced. - noise_map_2d - The 2D noise-map array whose values are replaced. - target_signal_to_noise - The target signal-to-noise the noise-map valueus are changed to. - - Returns - ------- - ndarray - The 2D noise-map with values changed. - - Examples - -------- - image_2d = np.ones((5,5)) - image_2d[2,2] = -1.0 - noise_map_2d = np.ones((5,5)) - - noise_map_2d_replaced = replace_noise_map_2d_values_where_image_2d_values_are_negative( - image_2d=image_2d, noise_map_2d=noise_map_2d, target_signal_to_noise=2.0): - """ - for y in range(image_2d.shape[0]): - for x in range(image_2d.shape[1]): - if image_2d[y, x] < 0.0: - absolute_signal_to_noise = np.abs(image_2d[y, x]) / noise_map_2d[y, x] - if absolute_signal_to_noise >= target_signal_to_noise: - noise_map_2d[y, x] = np.abs(image_2d[y, x]) / target_signal_to_noise - - return noise_map_2d - - @numba_util.jit() def index_2d_for_index_slim_from(indexes_slim: np.ndarray, shape_native) -> np.ndarray: """ diff --git a/test_autoarray/structures/arrays/test_array_2d_util.py b/test_autoarray/structures/arrays/test_array_2d_util.py index 0fe0d0a01..91a3fc70e 100644 --- a/test_autoarray/structures/arrays/test_array_2d_util.py +++ b/test_autoarray/structures/arrays/test_array_2d_util.py @@ -286,95 +286,6 @@ def test__resized_array_2d_from__padding_with_new_origin(): ).all() -def test__replace_noise_map_2d_values_where_image_2d_values_are_negative(): - image_2d = np.ones(shape=(2, 2)) - - noise_map_2d = np.array([[1.0, 2.0], [3.0, 4.0]]) - - noise_map_2d = ( - util.array_2d.replace_noise_map_2d_values_where_image_2d_values_are_negative( - image_2d=image_2d, noise_map_2d=noise_map_2d, target_signal_to_noise=1.0 - ) - ) - - assert (noise_map_2d == noise_map_2d).all() - - image_2d = -1.0 * np.ones(shape=(2, 2)) - - noise_map_2d = np.array([[1.0, 0.5], [0.25, 0.125]]) - - noise_map_2d = ( - util.array_2d.replace_noise_map_2d_values_where_image_2d_values_are_negative( - image_2d=image_2d, noise_map_2d=noise_map_2d, target_signal_to_noise=10.0 - ) - ) - - assert (noise_map_2d == noise_map_2d).all() - - noise_map_2d = ( - util.array_2d.replace_noise_map_2d_values_where_image_2d_values_are_negative( - image_2d=image_2d, noise_map_2d=noise_map_2d, target_signal_to_noise=4.0 - ) - ) - - assert (noise_map_2d == np.array([[1.0, 0.5], [0.25, 0.25]])).all() - - noise_map_2d = np.array([[1.0, 0.5], [0.25, 0.125]]) - - noise_map_2d = ( - util.array_2d.replace_noise_map_2d_values_where_image_2d_values_are_negative( - image_2d=image_2d, noise_map_2d=noise_map_2d, target_signal_to_noise=2.0 - ) - ) - - assert (noise_map_2d == np.array([[1.0, 0.5], [0.5, 0.5]])).all() - - noise_map_2d = np.array([[1.0, 0.5], [0.25, 0.125]]) - - noise_map_2d = ( - util.array_2d.replace_noise_map_2d_values_where_image_2d_values_are_negative( - image_2d=image_2d, noise_map_2d=noise_map_2d, target_signal_to_noise=1.0 - ) - ) - - assert (noise_map_2d == np.array([[1.0, 1.0], [1.0, 1.0]])).all() - - noise_map_2d = np.array([[1.0, 0.5], [0.25, 0.125]]) - - noise_map_2d = ( - util.array_2d.replace_noise_map_2d_values_where_image_2d_values_are_negative( - image_2d=image_2d, noise_map_2d=noise_map_2d, target_signal_to_noise=0.5 - ) - ) - - assert (noise_map_2d == np.array([[2.0, 2.0], [2.0, 2.0]])).all() - - -def test__same_as_above__image_not_all_negative(): - image_2d = np.array([[1.0, -2.0], [5.0, -4.0]]) - - noise_map_2d = np.array([[3.0, 1.0], [4.0, 8.0]]) - - noise_map_2d = ( - util.array_2d.replace_noise_map_2d_values_where_image_2d_values_are_negative( - image_2d=image_2d, noise_map_2d=noise_map_2d, target_signal_to_noise=1.0 - ) - ) - - assert (noise_map_2d == np.array([[3.0, 2.0], [4.0, 8.0]])).all() - - image_2d = np.array([[-10.0, -20.0], [100.0, -30.0]]) - - noise_map_2d = np.array([[1.0, 2.0], [40.0, 3.0]]) - - noise_map_2d = ( - util.array_2d.replace_noise_map_2d_values_where_image_2d_values_are_negative( - image_2d=image_2d, noise_map_2d=noise_map_2d, target_signal_to_noise=5.0 - ) - ) - - assert (noise_map_2d == np.array([[2.0, 4.0], [40.0, 6.0]])).all() - def test__index_2d_for_index_slim_from(): indexes_1d = np.array([0, 1, 2, 3, 4, 5, 6, 7, 8]) From 648a815ea51765005930245fff60b78ae45b320c Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Wed, 18 Jun 2025 20:44:11 +0100 Subject: [PATCH 03/19] index_2d_for_index_slim_from no longer uses numba --- autoarray/structures/arrays/array_2d_util.py | 11 +++++++---- 1 file changed, 7 insertions(+), 4 deletions(-) diff --git a/autoarray/structures/arrays/array_2d_util.py b/autoarray/structures/arrays/array_2d_util.py index a9be13fe6..06271c755 100644 --- a/autoarray/structures/arrays/array_2d_util.py +++ b/autoarray/structures/arrays/array_2d_util.py @@ -376,11 +376,14 @@ def index_2d_for_index_slim_from(indexes_slim: np.ndarray, shape_native) -> np.n indexes_slim = np.array([0, 1, 2, 5]) indexes_2d = index_2d_for_index_slim_from(indexes_slim=indexes_slim, shape=(3,3)) """ - index_2d_for_index_slim = np.zeros((indexes_slim.shape[0], 2)) + # Calculate row indices by integer division by number of columns + rows = indexes_slim // shape_native[1] - for i, index_slim in enumerate(indexes_slim): - index_2d_for_index_slim[i, 0] = int(index_slim / shape_native[1]) - index_2d_for_index_slim[i, 1] = int(index_slim % shape_native[1]) + # Calculate column indices by modulo number of columns + cols = indexes_slim % shape_native[1] + + # Stack rows and cols horizontally into shape (N, 2) + index_2d_for_index_slim = np.vstack((rows, cols)).T return index_2d_for_index_slim From 4a8814a0a5e0c7819ec420c34b01713705e7bd43 Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Wed, 18 Jun 2025 20:45:10 +0100 Subject: [PATCH 04/19] index_slim_for_index_2d_from does not use numba --- autoarray/structures/arrays/array_2d_util.py | 12 ++---------- 1 file changed, 2 insertions(+), 10 deletions(-) diff --git a/autoarray/structures/arrays/array_2d_util.py b/autoarray/structures/arrays/array_2d_util.py index 06271c755..8cf1f0cf2 100644 --- a/autoarray/structures/arrays/array_2d_util.py +++ b/autoarray/structures/arrays/array_2d_util.py @@ -342,8 +342,6 @@ def resized_array_2d_from( return resized_array - -@numba_util.jit() def index_2d_for_index_slim_from(indexes_slim: np.ndarray, shape_native) -> np.ndarray: """ For pixels on a native 2D array of shape (total_y_pixels, total_x_pixels), this array maps the slimmed 1D pixel @@ -387,8 +385,6 @@ def index_2d_for_index_slim_from(indexes_slim: np.ndarray, shape_native) -> np.n return index_2d_for_index_slim - -@numba_util.jit() def index_slim_for_index_2d_from(indexes_2d: np.ndarray, shape_native) -> np.ndarray: """ For pixels on a native 2D array of shape (total_y_pixels, total_x_pixels), this array maps the 2D pixel indexes to @@ -421,12 +417,8 @@ def index_slim_for_index_2d_from(indexes_2d: np.ndarray, shape_native) -> np.nda indexes_2d = np.array([[0,0], [1,0], [2,0], [2,2]]) indexes_flat = index_flat_for_index_2d_from(indexes_2d=indexes_2d, shape=(3,3)) """ - index_slim_for_index_native_2d = np.zeros(indexes_2d.shape[0]) - - for i in range(indexes_2d.shape[0]): - index_slim_for_index_native_2d[i] = int( - (indexes_2d[i, 0]) * shape_native[1] + indexes_2d[i, 1] - ) + # Calculate 1D indexes as row_index * number_of_columns + col_index + index_slim_for_index_native_2d = indexes_2d[:, 0] * shape_native[1] + indexes_2d[:, 1] return index_slim_for_index_native_2d From d69e39f8e9dd5e406fffc2b5b2664c076df3086d Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Wed, 18 Jun 2025 20:46:38 +0100 Subject: [PATCH 05/19] remove array_2d_slim_complex_from --- autoarray/structures/arrays/array_2d_util.py | 45 ------------------- .../structures/arrays/test_array_2d_util.py | 18 -------- 2 files changed, 63 deletions(-) diff --git a/autoarray/structures/arrays/array_2d_util.py b/autoarray/structures/arrays/array_2d_util.py index 8cf1f0cf2..915820f07 100644 --- a/autoarray/structures/arrays/array_2d_util.py +++ b/autoarray/structures/arrays/array_2d_util.py @@ -558,51 +558,6 @@ def array_2d_via_indexes_from( ) -@numba_util.jit() -def array_2d_slim_complex_from( - array_2d_native: np.ndarray, - mask: np.ndarray, -) -> np.ndarray: - """ - For a 2D array and mask, map the values of all unmasked pixels to a 1D array. - - The pixel coordinate origin is at the top left corner of the 2D array and goes right-wards and downwards. - - For example, for an array of shape (3,3) and where all pixels are unmasked: - - - pixel [0,0] of the 2D array will correspond to index 0 of the 1D array. - - pixel [0,1] of the 2D array will correspond to index 1 of the 1D array. - - pixel [1,0] of the 2D array will correspond to index 3 of the 1D array. - - pixel [2,0] of the 2D array will correspond to index 6 of the 1D array. - - Parameters - ---------- - array_2d_native - A 2D array of values on the dimensions of the grid. - mask - A 2D array of bools, where `False` values mean unmasked and are included in the mapping. - array_2d - The 2D array of values which are mapped to a 1D array. - - Returns - ------- - ndarray - A 1D array of values mapped from the 2D array with dimensions (total_unmasked_pixels). - """ - - total_pixels = np.sum(~mask) - - array_1d = 0 + 0j * np.zeros(shape=total_pixels) - index = 0 - - for y in range(mask.shape[0]): - for x in range(mask.shape[1]): - if not mask[y, x]: - array_1d[index] = array_2d_native[y, x] - index += 1 - - return array_1d - @numba_util.jit() def array_2d_native_complex_via_indexes_from( diff --git a/test_autoarray/structures/arrays/test_array_2d_util.py b/test_autoarray/structures/arrays/test_array_2d_util.py index 91a3fc70e..144c8d39d 100644 --- a/test_autoarray/structures/arrays/test_array_2d_util.py +++ b/test_autoarray/structures/arrays/test_array_2d_util.py @@ -431,24 +431,6 @@ def test__array_2d_slim_from(): assert (array_2d_slim == np.array([2, 4, 5, 6, 8])).all() -def test__array_2d_slim_from__complex_array(): - array_2d = np.array( - [ - [1 + 1j, 2 + 2j, 3 + 3], - [4 + 4j, 5 + 5j, 6 + 6j], - [7 + 7j, 8 + 8j, 9 + 9j], - ] - ) - - mask = np.array([[True, True, True], [True, False, True], [True, True, True]]) - - array_2d_slim = util.array_2d.array_2d_slim_complex_from( - mask=mask, - array_2d_native=array_2d, - ) - - assert (array_2d_slim == np.array([5 + 5j])).all() - def test__array_2d_native_from(): array_2d_slim = np.array([1.0, 2.0, 3.0, 4.0]) From 3832bca6a6b3e43c073d6489d7e2f4d9446841d2 Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Wed, 18 Jun 2025 20:46:43 +0100 Subject: [PATCH 06/19] remove array_2d_native_complex_via_indexes_from --- autoarray/structures/arrays/array_2d_util.py | 18 ------------------ .../structures/arrays/test_array_2d_util.py | 15 --------------- 2 files changed, 33 deletions(-) diff --git a/autoarray/structures/arrays/array_2d_util.py b/autoarray/structures/arrays/array_2d_util.py index 915820f07..32b33a79b 100644 --- a/autoarray/structures/arrays/array_2d_util.py +++ b/autoarray/structures/arrays/array_2d_util.py @@ -556,21 +556,3 @@ def array_2d_via_indexes_from( return ( jnp.zeros(shape).at[tuple(native_index_for_slim_index_2d.T)].set(array_2d_slim) ) - - - -@numba_util.jit() -def array_2d_native_complex_via_indexes_from( - array_2d_slim: np.ndarray, - shape_native: Tuple[int, int], - native_index_for_slim_index_2d: np.ndarray, -) -> np.ndarray: - array_2d = 0 + 0j * np.zeros(shape_native) - - for slim_index in range(len(native_index_for_slim_index_2d)): - array_2d[ - native_index_for_slim_index_2d[slim_index, 0], - native_index_for_slim_index_2d[slim_index, 1], - ] = array_2d_slim[slim_index] - - return array_2d diff --git a/test_autoarray/structures/arrays/test_array_2d_util.py b/test_autoarray/structures/arrays/test_array_2d_util.py index 144c8d39d..8b8056f3b 100644 --- a/test_autoarray/structures/arrays/test_array_2d_util.py +++ b/test_autoarray/structures/arrays/test_array_2d_util.py @@ -477,18 +477,3 @@ def test__array_2d_native_from(): ) ).all() - -def test__array_2d_native_from__compelx_array(): - array_2d_slim = np.array( - [1.0 + 1j, 2.0 + 2j, 3.0 + 3j, 4.0 + 4j], dtype="complex128" - ) - - array_2d = util.array_2d.array_2d_native_complex_via_indexes_from( - array_2d_slim=array_2d_slim, - shape_native=(2, 2), - native_index_for_slim_index_2d=np.array( - [[0, 0], [0, 1], [1, 0], [1, 1]], dtype="int" - ), - ) - - assert (array_2d == np.array([[1.0 + 1j, 2.0 + 2j], [3.0 + 3j, 4.0 + 4j]])).all() From 783ff848b1af41d767afd4ed749ef03065e12bc1 Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Wed, 18 Jun 2025 20:51:22 +0100 Subject: [PATCH 07/19] remove numba from two more functions ing rid_2d_util --- autoarray/structures/grids/grid_2d_util.py | 11 +++++------ 1 file changed, 5 insertions(+), 6 deletions(-) diff --git a/autoarray/structures/grids/grid_2d_util.py b/autoarray/structures/grids/grid_2d_util.py index b547fcb94..f8f78e0f0 100644 --- a/autoarray/structures/grids/grid_2d_util.py +++ b/autoarray/structures/grids/grid_2d_util.py @@ -395,8 +395,6 @@ def grid_2d_via_shape_native_from( origin=origin, ) - -@numba_util.jit() def _radial_projected_shape_slim_from( extent: np.ndarray, centre: Tuple[float, float], @@ -471,7 +469,6 @@ def _radial_projected_shape_slim_from( return int((scaled_distance / pixel_scale)) + 1 -@numba_util.jit() def grid_scaled_2d_slim_radial_projected_from( extent: np.ndarray, centre: Tuple[float, float], @@ -562,9 +559,11 @@ def grid_scaled_2d_slim_radial_projected_from( radii = centre[1] - for slim_index in range(shape_slim): - grid_scaled_2d_slim_radii[slim_index, 1] = radii - radii += pixel_scale + # Create an array of radii values spaced by pixel_scale + radii_array = radii + pixel_scale * np.arange(shape_slim) + + # Assign all values at once to the second column (index 1) + grid_scaled_2d_slim_radii[:, 1] = radii_array return grid_scaled_2d_slim_radii + 1e-6 From ab45cf50e84fb67ca52e9e9879997706f5d375a4 Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Wed, 18 Jun 2025 20:52:25 +0100 Subject: [PATCH 08/19] remove grid_2d_slim_upscaled_from --- autoarray/structures/grids/grid_2d_util.py | 52 --------- .../structures/grids/test_grid_2d_util.py | 106 +----------------- 2 files changed, 1 insertion(+), 157 deletions(-) diff --git a/autoarray/structures/grids/grid_2d_util.py b/autoarray/structures/grids/grid_2d_util.py index f8f78e0f0..00c67cd7d 100644 --- a/autoarray/structures/grids/grid_2d_util.py +++ b/autoarray/structures/grids/grid_2d_util.py @@ -738,58 +738,6 @@ def grid_2d_native_from( return jnp.stack((grid_2d_native_y, grid_2d_native_x), axis=-1) -@numba_util.jit() -def grid_2d_slim_upscaled_from( - grid_slim: np.ndarray, upscale_factor: int, pixel_scales: ty.PixelScales -) -> np.ndarray: - """ - From an input slimmed 2D grid, return an upscaled slimmed 2D grid where (y,x) coordinates are added at an - upscaled resolution to each grid coordinate. - - Parameters - ---------- - grid_slim - The slimmed grid of (y,x) coordinates over which a square uniform grid is overlaid. - upscale_factor - The upscaled resolution at which the new grid coordinates are computed. - pixel_scales - The pixel scale of the uniform grid that laid over the irregular grid of (y,x) coordinates. - """ - - grid_2d_slim_upscaled = np.zeros(shape=(grid_slim.shape[0] * upscale_factor**2, 2)) - - upscale_index = 0 - - y_upscale_half = pixel_scales[0] / 2 - y_upscale_step = pixel_scales[0] / upscale_factor - - x_upscale_half = pixel_scales[1] / 2 - x_upscale_step = pixel_scales[1] / upscale_factor - - for slim_index in range(grid_slim.shape[0]): - y_grid = grid_slim[slim_index, 0] - x_grid = grid_slim[slim_index, 1] - - for y in range(upscale_factor): - for x in range(upscale_factor): - grid_2d_slim_upscaled[upscale_index, 0] = ( - y_grid - + y_upscale_half - - y * y_upscale_step - - (y_upscale_step / 2.0) - ) - grid_2d_slim_upscaled[upscale_index, 1] = ( - x_grid - - x_upscale_half - + x * x_upscale_step - + (x_upscale_step / 2.0) - ) - - upscale_index += 1 - - return grid_2d_slim_upscaled - - def grid_2d_of_points_within_radius( radius: float, centre: Tuple[float, float], grid_2d: np.ndarray ): diff --git a/test_autoarray/structures/grids/test_grid_2d_util.py b/test_autoarray/structures/grids/test_grid_2d_util.py index 463aa4e7b..4adc6c020 100644 --- a/test_autoarray/structures/grids/test_grid_2d_util.py +++ b/test_autoarray/structures/grids/test_grid_2d_util.py @@ -460,108 +460,4 @@ def test__grid_2d_native_from(): [[-1.0, -1.0], [-2.0, -2.0], [0.0, 0.0], [-3.0, -3.0]], ] ) - ).all() - - grid_slim = np.array( - [ - [1.0, 1.0], - [1.0, 1.0], - [1.0, 1.0], - [1.0, 1.0], - [2.0, 2.0], - [2.0, 2.0], - [2.0, 2.0], - [2.0, 2.0], - [3.0, 3.0], - [3.0, 3.0], - [3.0, 3.0], - [4.0, 4.0], - ] - ) - - -def test__grid_2d_slim_upscaled_from(): - grid_slim = np.array([[1.0, 1.0]]) - - grid_upscaled_2d = aa.util.grid_2d.grid_2d_slim_upscaled_from( - grid_slim=grid_slim, upscale_factor=1, pixel_scales=(2.0, 2.0) - ) - - assert (grid_upscaled_2d == np.array([[1.0, 1.0]])).all() - - grid_upscaled_2d = aa.util.grid_2d.grid_2d_slim_upscaled_from( - grid_slim=grid_slim, upscale_factor=2, pixel_scales=(2.0, 2.0) - ) - - assert ( - grid_upscaled_2d == np.array([[1.5, 0.5], [1.5, 1.5], [0.5, 0.5], [0.5, 1.5]]) - ).all() - - grid_slim = np.array([[1.0, 1.0], [1.0, 3.0]]) - - grid_upscaled_2d = aa.util.grid_2d.grid_2d_slim_upscaled_from( - grid_slim=grid_slim, upscale_factor=2, pixel_scales=(2.0, 2.0) - ) - - assert ( - grid_upscaled_2d - == np.array( - [ - [1.5, 0.5], - [1.5, 1.5], - [0.5, 0.5], - [0.5, 1.5], - [1.5, 2.5], - [1.5, 3.5], - [0.5, 2.5], - [0.5, 3.5], - ] - ) - ).all() - - grid_slim = np.array([[1.0, 1.0], [3.0, 1.0]]) - - grid_upscaled_2d = aa.util.grid_2d.grid_2d_slim_upscaled_from( - grid_slim=grid_slim, upscale_factor=2, pixel_scales=(2.0, 2.0) - ) - - assert ( - grid_upscaled_2d - == np.array( - [ - [1.5, 0.5], - [1.5, 1.5], - [0.5, 0.5], - [0.5, 1.5], - [3.5, 0.5], - [3.5, 1.5], - [2.5, 0.5], - [2.5, 1.5], - ] - ) - ).all() - - grid_slim = np.array([[1.0, 1.0]]) - - grid_upscaled_2d = aa.util.grid_2d.grid_2d_slim_upscaled_from( - grid_slim=grid_slim, upscale_factor=2, pixel_scales=(3.0, 2.0) - ) - - assert ( - grid_upscaled_2d - == np.array([[1.75, 0.5], [1.75, 1.5], [0.25, 0.5], [0.25, 1.5]]) - ).all() - - grid_upscaled_2d = aa.util.grid_2d.grid_2d_slim_upscaled_from( - grid_slim=grid_slim, upscale_factor=3, pixel_scales=(2.0, 2.0) - ) - - assert grid_upscaled_2d[0] == pytest.approx(np.array([1.666, 0.333]), 1.0e-2) - assert grid_upscaled_2d[1] == pytest.approx(np.array([1.666, 1.0]), 1.0e-2) - assert grid_upscaled_2d[2] == pytest.approx(np.array([1.666, 1.666]), 1.0e-2) - assert grid_upscaled_2d[3] == pytest.approx(np.array([1.0, 0.333]), 1.0e-2) - assert grid_upscaled_2d[4] == pytest.approx(np.array([1.0, 1.0]), 1.0e-2) - assert grid_upscaled_2d[5] == pytest.approx(np.array([1.0, 1.666]), 1.0e-2) - assert grid_upscaled_2d[6] == pytest.approx(np.array([0.333, 0.333]), 1.0e-2) - assert grid_upscaled_2d[7] == pytest.approx(np.array([0.333, 1.0]), 1.0e-2) - assert grid_upscaled_2d[8] == pytest.approx(np.array([0.333, 1.666]), 1.0e-2) + ).all() \ No newline at end of file From 8bb449a1fd18ac75e6ae3a971f9c0bdba4405827 Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Wed, 18 Jun 2025 20:59:44 +0100 Subject: [PATCH 09/19] remove native_sub_index_for_slim_sub_index_2d_from --- .../over_sampling/over_sample_util.py | 83 ---------- .../operators/over_sampling/over_sampler.py | 59 +------ autoarray/structures/arrays/array_2d_util.py | 1 - autoarray/util/cholesky_funcs.py | 100 ----------- autoarray/util/fnnls.py | 155 ------------------ 5 files changed, 1 insertion(+), 397 deletions(-) delete mode 100644 autoarray/util/cholesky_funcs.py delete mode 100644 autoarray/util/fnnls.py diff --git a/autoarray/operators/over_sampling/over_sample_util.py b/autoarray/operators/over_sampling/over_sample_util.py index d3488a457..7f7f11312 100644 --- a/autoarray/operators/over_sampling/over_sample_util.py +++ b/autoarray/operators/over_sampling/over_sample_util.py @@ -78,89 +78,6 @@ def total_sub_pixels_2d_from(sub_size: np.ndarray) -> int: return int(np.sum(sub_size**2)) -@numba_util.jit() -def native_sub_index_for_slim_sub_index_2d_from( - mask_2d: np.ndarray, sub_size: np.ndarray -) -> np.ndarray: - """ - Returns an array of shape [total_unmasked_pixels*sub_size] that maps every unmasked sub-pixel to its - corresponding native 2D pixel using its (y,x) pixel indexes. - - For example, for the following ``Mask2D`` for ``sub_size=1``: - - :: - [[True, True, True, True] - [True, False, False, True], - [True, False, True, True], - [True, True, True, True]] - - This has three unmasked (``False`` values) which have the ``slim`` indexes: - - :: - [0, 1, 2] - - The array ``native_index_for_slim_index_2d`` is therefore: - - :: - [[1,1], [1,2], [2,1]] - - For a ``Mask2D`` with ``sub_size=2`` each unmasked ``False`` entry is split into a sub-pixel of size 2x2 and - there are therefore 12 ``slim`` indexes: - - :: - [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11] - - The array ``native_index_for_slim_index_2d`` is therefore: - - :: - [[2,2], [2,3], [2,4], [2,5], [3,2], [3,3], [3,4], [3,5], [4,2], [4,3], [5,2], [5,3]] - - Parameters - ---------- - mask_2d - A 2D array of bools, where `False` values are unmasked. - sub_size - The size of the sub-grid in each mask pixel. - - Returns - ------- - ndarray - An array that maps pixels from a slimmed array of shape [total_unmasked_pixels*sub_size] to its native array - of shape [total_pixels*sub_size, total_pixels*sub_size]. - - Examples - -------- - mask_2d = np.array([[True, True, True], - [True, False, True] - [True, True, True]]) - - sub_native_index_for_sub_slim_index_2d = sub_native_index_for_sub_slim_index_via_mask_2d_from(mask_2d=mask_2d, sub_size=1) - """ - - total_sub_pixels = total_sub_pixels_2d_from(sub_size=sub_size) - sub_native_index_for_sub_slim_index_2d = np.zeros(shape=(total_sub_pixels, 2)) - - slim_index = 0 - sub_slim_index = 0 - - for y in range(mask_2d.shape[0]): - for x in range(mask_2d.shape[1]): - if not mask_2d[y, x]: - sub = sub_size[slim_index] - - for y1 in range(sub): - for x1 in range(sub): - sub_native_index_for_sub_slim_index_2d[sub_slim_index, :] = ( - (y * sub) + y1, - (x * sub) + x1, - ) - sub_slim_index += 1 - - slim_index += 1 - - return sub_native_index_for_sub_slim_index_2d - - @numba_util.jit() def slim_index_for_sub_slim_index_via_mask_2d_from( mask_2d: np.ndarray, sub_size: np.ndarray diff --git a/autoarray/operators/over_sampling/over_sampler.py b/autoarray/operators/over_sampling/over_sampler.py index 5d9a99871..268cc7be9 100644 --- a/autoarray/operators/over_sampling/over_sampler.py +++ b/autoarray/operators/over_sampling/over_sampler.py @@ -128,7 +128,7 @@ def __init__(self, mask: Mask2D, sub_size: Union[int, Array2D]): based on the sub-grid sizes. The over sampling class has functions dedicated to mapping between the sub-grid and pixel-grid, for example - `sub_mask_native_for_sub_mask_slim` and `slim_for_sub_slim`. + `slim_for_sub_slim`. The class `OverSampling` is used for the high level API, whereby this is where users input their preferred over-sampling configuration. This class, `OverSampler`, contains the functionality @@ -258,63 +258,6 @@ def binned_array_2d_from(self, array: Array2D) -> "Array2D": mask=self.mask, ) - @cached_property - def sub_mask_native_for_sub_mask_slim(self) -> np.ndarray: - """ - Derives a 1D ``ndarray`` which maps every subgridded 1D ``slim`` index of the ``Mask2D`` to its - subgridded 2D ``native`` index. - - For example, for the following ``Mask2D`` for ``sub_size=1``: - - :: - [[True, True, True, True] - [True, False, False, True], - [True, False, True, True], - [True, True, True, True]] - - This has three unmasked (``False`` values) which have the ``slim`` indexes: - - :: - [0, 1, 2] - - The array ``sub_mask_native_for_sub_mask_slim`` is therefore: - - :: - [[1,1], [1,2], [2,1]] - - For a ``Mask2D`` with ``sub_size=2`` each unmasked ``False`` entry is split into a sub-pixel of size 2x2 and - there are therefore 12 ``slim`` indexes: - - :: - [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11] - - The array ``native_for_slim`` is therefore: - - :: - [[2,2], [2,3], [2,4], [2,5], [3,2], [3,3], [3,4], [3,5], [4,2], [4,3], [5,2], [5,3]] - - Examples - -------- - - .. code-block:: python - - import autoarray as aa - - mask_2d = aa.Mask2D( - mask=[[True, True, True, True] - [True, False, False, True], - [True, False, True, True], - [True, True, True, True]] - pixel_scales=1.0, - ) - - derive_indexes_2d = aa.DeriveIndexes2D(mask=mask_2d) - - print(derive_indexes_2d.sub_mask_native_for_sub_mask_slim) - """ - return over_sample_util.native_sub_index_for_slim_sub_index_2d_from( - mask_2d=self.mask.array, sub_size=self.sub_size.array - ).astype("int") @cached_property def slim_for_sub_slim(self) -> np.ndarray: diff --git a/autoarray/structures/arrays/array_2d_util.py b/autoarray/structures/arrays/array_2d_util.py index 32b33a79b..92e99d43a 100644 --- a/autoarray/structures/arrays/array_2d_util.py +++ b/autoarray/structures/arrays/array_2d_util.py @@ -6,7 +6,6 @@ if TYPE_CHECKING: from autoarray.mask.mask_2d import Mask2D -from autoarray import numba_util from autoarray.mask import mask_2d_util from autoarray import exc diff --git a/autoarray/util/cholesky_funcs.py b/autoarray/util/cholesky_funcs.py deleted file mode 100644 index bd211eeb5..000000000 --- a/autoarray/util/cholesky_funcs.py +++ /dev/null @@ -1,100 +0,0 @@ -import numpy as np -from scipy import linalg -import math -import time -from autoarray import numba_util - - -@numba_util.jit() -def _choldowndate(U, x): - n = x.size - for k in range(n - 1): - Ukk = U[k, k] - xk = x[k] - r = math.sqrt(Ukk**2 - xk**2) - c = r / Ukk - s = xk / Ukk - U[k, k] = r - U[k, k + 1 :] = (U[k, (k + 1) :] - s * x[k + 1 :]) / c - x[k + 1 :] = c * x[k + 1 :] - s * U[k, k + 1 :] - - k = n - 1 - U[k, k] = math.sqrt(U[k, k] ** 2 - x[k] ** 2) - return U - - -@numba_util.jit() -def _cholupdate(U, x): - n = x.size - for k in range(n - 1): - Ukk = U[k, k] - xk = x[k] - - r = np.sqrt(Ukk**2 + xk**2) - - c = r / Ukk - s = xk / Ukk - U[k, k] = r - - U[k, k + 1 :] = (U[k, (k + 1) :] + s * x[k + 1 :]) / c - x[k + 1 :] = c * x[k + 1 :] - s * U[k, k + 1 :] - - k = n - 1 - U[k, k] = np.sqrt(U[k, k] ** 2 + x[k] ** 2) - - return U - - -def cholinsert(U, index, x): - S = np.insert(np.insert(U, index, 0, axis=0), index, 0, axis=1) - - S[:index, index] = S12 = linalg.solve_triangular( - U[:index, :index], x[:index], trans=1, lower=False, overwrite_b=True - ) - - S[index, index] = s22 = math.sqrt(x[index] - S12.dot(S12)) - - if index == U.shape[0]: - return S - else: - S[index, index + 1 :] = S23 = (x[index + 1 :] - S12.T @ U[:index, index:]) / s22 - _choldowndate(S[index + 1 :, index + 1 :], S23) # S33 - return S - - -def cholinsertlast(U, x): - """ - Update the Cholesky matrix U by inserting a vector at the end of the matrix - Inserting a vector to the end of U doesn't require _cholupdate, so save some time. - It's a special case of `cholinsert` (as shown above, if index == U.shape[0]) - As in current Cholesky scheme implemented in fnnls, we only use this kind of insertion, so I - separate it out from the `cholinsert`. - """ - index = U.shape[0] - - S = np.insert(np.insert(U, index, 0, axis=0), index, 0, axis=1) - - S[:index, index] = S12 = linalg.solve_triangular( - U[:index, :index], x[:index], trans=1, lower=False, overwrite_b=True - ) - - S[index, index] = s22 = math.sqrt(x[index] - S12.dot(S12)) - - return S - - -def choldeleteindexes(U, indexes): - indexes = sorted(indexes, reverse=True) - - for index in indexes: - L = np.delete(np.delete(U, index, axis=0), index, axis=1) - - # If the deleted index is at the end of matrix, then we do not need to update the U. - - if index == L.shape[0]: - U = L - else: - _cholupdate(L[index:, index:], U[index, index + 1 :]) - U = L - - return U diff --git a/autoarray/util/fnnls.py b/autoarray/util/fnnls.py deleted file mode 100644 index 3f49c1f2d..000000000 --- a/autoarray/util/fnnls.py +++ /dev/null @@ -1,155 +0,0 @@ -import numpy as np -from scipy import linalg as slg - -from autoarray.util.cholesky_funcs import cholinsertlast, choldeleteindexes - -from autoarray import exc - -""" - This file contains functions use the Bro & Jong (1997) algorithm to solve the non-negative least - square problem. The `fnnls and fix_constraint` is orginally copied from - "https://github.com/jvendrow/fnnls". - For our purpose in PyAutoArray, we create `fnnls_modefied` to take ZTZ and ZTx as inputs directly. - Furthermore, we add two functions `fnnls_Cholesky and fix_constraint_Cholesky` to realize a scheme - that solves the lstsq problem in the algorithm by Cholesky factorisation. For ~ 1000 free - parameters, we see a speed up by 2 times and should be more for more parameters. - We have also noticed that by setting the P_initial to be `sla.solve(ZTZ, ZTx, assume_a='pos') > 0` - will speed up our task (~ 1000 free parameters) by ~ 3 times as it significantly reduces the - iteration time. -""" - - -def fnnls_cholesky( - ZTZ, - ZTx, - P_initial=np.zeros(0, dtype=int), -): - """ - Similar to fnnls, but use solving the lstsq problem by updating Cholesky factorisation. - """ - - lstsq = lambda A, x: slg.solve( - A, - x, - assume_a="pos", - overwrite_a=True, - overwrite_b=True, - ) - - n = np.shape(ZTZ)[0] - epsilon = 2.2204e-16 - tolerance = epsilon * n - max_repetitions = 3 - no_update = 0 - loop_count = 0 - loop_count2 = 0 - - P = np.zeros(n, dtype=bool) - P[P_initial] = True - d = np.zeros(n) - w = ZTx - (ZTZ) @ d - s_chol = np.zeros(n) - - if P_initial.shape[0] != 0: - P_number = np.arange(len(P), dtype="int") - P_inorder = P_number[P_initial] - s_chol[P] = lstsq((ZTZ)[P][:, P], (ZTx)[P]) - d = s_chol.clip(min=0) - else: - P_inorder = np.array([], dtype="int") - - # P_inorder is similar as P. They are both used to select solutions in the passive set. - # P_inorder saves the `indexes` of those passive solutions. - # P saves [True/False] for all solutions. True indicates a solution in the passive set while False - # indicates it's in the active set. - # The benifit of P_inorder is that we are able to not only select out solutions in the passive set - # and can sort them in the order of added to the passive set. This will make updating the - # Cholesky factorisation simpler and thus save time. - - while (not np.all(P)) and np.max(w[~P]) > tolerance: - # make copy of passive set to check for change at end of loop - - current_P = P.copy() - idmax = np.argmax(w * ~P) - P_inorder = np.append(P_inorder, int(idmax)) - - if loop_count == 0: - # We need to initialize the Cholesky factorisation, U, for the first loop. - U = slg.cholesky(ZTZ[P_inorder][:, P_inorder]) - else: - U = cholinsertlast(U, ZTZ[idmax][P_inorder]) - - # solve the lstsq problem by cho_solve - - s_chol[P_inorder] = slg.cho_solve((U, False), ZTx[P_inorder]) - - P[idmax] = True - while np.any(P) and np.min(s_chol[P]) <= tolerance: - s_chol, d, P, P_inorder, U = fix_constraint_cholesky( - ZTx=ZTx, - s_chol=s_chol, - d=d, - P=P, - P_inorder=P_inorder, - U=U, - tolerance=tolerance, - ) - - loop_count2 += 1 - if loop_count2 > 10000: - raise RuntimeError - - d = s_chol.copy() - w = ZTx - (ZTZ) @ d - loop_count += 1 - - if loop_count > 10000: - raise RuntimeError - - if np.all(current_P == P): - no_update += 1 - else: - no_update = 0 - - if no_update >= max_repetitions: - break - - return d - - -def fix_constraint_cholesky(ZTx, s_chol, d, P, P_inorder, U, tolerance): - """ - Similar to fix_constraint, but solve the lstsq by Cholesky factorisation. - If this function is called, it means some solutions in the current passive sets needed to be - taken out and put into the active set. - So, this function involves 3 procedure: - 1. Identifying what solutions should be taken out of the current passive set. - 2. Updating the P, P_inorder and the Cholesky factorisation U. - 3. Solving the lstsq by using the new Cholesky factorisation U. - As some solutions are taken out from the passive set, the Cholesky factorisation needs to be - updated by choldeleteindexes. To realize that, we call the `choldeleteindexes` from - cholesky_funcs. - """ - q = P * (s_chol <= tolerance) - alpha = np.min(d[q] / (d[q] - s_chol[q])) - - # set d as close to s as possible while maintaining non-negativity - d = d + alpha * (s_chol - d) - - id_delete = np.where(d[P_inorder] <= tolerance)[0] - - U = choldeleteindexes(U, id_delete) # update the Cholesky factorisation - - P_inorder = np.delete(P_inorder, id_delete) # update the P_inorder - - P[d <= tolerance] = False # update the P - - # solve the lstsq problem by cho_solve - - if len(P_inorder): - # there could be a case where P_inorder is empty. - s_chol[P_inorder] = slg.cho_solve((U, False), ZTx[P_inorder]) - - s_chol[~P] = 0.0 # set solutions taken out of the passive set to be 0 - - return s_chol, d, P, P_inorder, U From 0f03df7730a1ee01b1ed3772e3ade9dc05e29d6d Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Wed, 18 Jun 2025 21:01:19 +0100 Subject: [PATCH 10/19] remove redudnant tests --- .../over_sampling/over_sample_util.py | 1 - .../over_sample/test_over_sample_util.py | 71 ------------------- .../over_sample/test_over_sampler.py | 20 ------ 3 files changed, 92 deletions(-) diff --git a/autoarray/operators/over_sampling/over_sample_util.py b/autoarray/operators/over_sampling/over_sample_util.py index 7f7f11312..d45aee825 100644 --- a/autoarray/operators/over_sampling/over_sample_util.py +++ b/autoarray/operators/over_sampling/over_sample_util.py @@ -78,7 +78,6 @@ def total_sub_pixels_2d_from(sub_size: np.ndarray) -> int: return int(np.sum(sub_size**2)) -@numba_util.jit() def slim_index_for_sub_slim_index_via_mask_2d_from( mask_2d: np.ndarray, sub_size: np.ndarray ) -> np.ndarray: diff --git a/test_autoarray/operators/over_sample/test_over_sample_util.py b/test_autoarray/operators/over_sample/test_over_sample_util.py index c0552d8a0..dc6237dc8 100644 --- a/test_autoarray/operators/over_sample/test_over_sample_util.py +++ b/test_autoarray/operators/over_sample/test_over_sample_util.py @@ -12,77 +12,6 @@ def test__total_sub_pixels_2d_from(): ) -def test__native_sub_index_for_slim_sub_index_2d_from(): - mask = np.array([[True, True, True], [True, False, True], [True, True, True]]) - - sub_mask_index_for_sub_mask_1d_index = ( - util.over_sample.native_sub_index_for_slim_sub_index_2d_from( - mask_2d=mask, sub_size=np.array([2]) - ) - ) - - assert ( - sub_mask_index_for_sub_mask_1d_index - == np.array([[2, 2], [2, 3], [3, 2], [3, 3]]) - ).all() - - mask = np.array([[True, False, True], [False, False, False], [True, False, True]]) - - sub_mask_index_for_sub_mask_1d_index = ( - util.over_sample.native_sub_index_for_slim_sub_index_2d_from( - mask_2d=mask, sub_size=np.array([2, 2, 2, 2, 2]) - ) - ) - - assert ( - sub_mask_index_for_sub_mask_1d_index - == np.array( - [ - [0, 2], - [0, 3], - [1, 2], - [1, 3], - [2, 0], - [2, 1], - [3, 0], - [3, 1], - [2, 2], - [2, 3], - [3, 2], - [3, 3], - [2, 4], - [2, 5], - [3, 4], - [3, 5], - [4, 2], - [4, 3], - [5, 2], - [5, 3], - ] - ) - ).all() - - mask = np.array( - [ - [True, True, True], - [True, False, True], - [True, True, True], - [True, True, False], - ] - ) - - sub_mask_index_for_sub_mask_1d_index = ( - util.over_sample.native_sub_index_for_slim_sub_index_2d_from( - mask_2d=mask, sub_size=np.array([2, 2]) - ) - ) - - assert ( - sub_mask_index_for_sub_mask_1d_index - == np.array([[2, 2], [2, 3], [3, 2], [3, 3], [6, 4], [6, 5], [7, 4], [7, 5]]) - ).all() - - def test__slim_index_for_sub_slim_index_via_mask_2d_from(): mask = np.array([[True, True, True], [True, False, True], [True, True, True]]) diff --git a/test_autoarray/operators/over_sample/test_over_sampler.py b/test_autoarray/operators/over_sample/test_over_sampler.py index a32b11e8f..e7d601614 100644 --- a/test_autoarray/operators/over_sample/test_over_sampler.py +++ b/test_autoarray/operators/over_sample/test_over_sampler.py @@ -91,26 +91,6 @@ def test__binned_array_2d_from(): assert binned_array_2d.slim == pytest.approx(np.array([1.0, 8.0]), 1.0e-4) -def test__sub_mask_index_for_sub_mask_1d_index(): - mask = aa.Mask2D( - mask=[[True, True, True], [True, False, False], [True, True, False]], - pixel_scales=1.0, - sub_size=2, - ) - - over_sampling = aa.OverSampler(mask=mask, sub_size=2) - - sub_mask_index_for_sub_mask_1d_index = ( - aa.util.over_sample.native_sub_index_for_slim_sub_index_2d_from( - mask_2d=np.array(mask), sub_size=np.array([2, 2, 2]) - ) - ) - - assert over_sampling.sub_mask_native_for_sub_mask_slim == pytest.approx( - sub_mask_index_for_sub_mask_1d_index, 1e-4 - ) - - def test__slim_index_for_sub_slim_index(): mask = aa.Mask2D( mask=[[True, False, True], [False, False, False], [True, False, False]], From 5a163af882ee9fec16308507b4a35fa2ac89c203 Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Wed, 18 Jun 2025 21:02:29 +0100 Subject: [PATCH 11/19] slim_index_for_sub_slim_index_via_mask_2d_from no longer uses numba --- .../over_sampling/over_sample_util.py | 23 +++++++------------ 1 file changed, 8 insertions(+), 15 deletions(-) diff --git a/autoarray/operators/over_sampling/over_sample_util.py b/autoarray/operators/over_sampling/over_sample_util.py index d45aee825..b86062d14 100644 --- a/autoarray/operators/over_sampling/over_sample_util.py +++ b/autoarray/operators/over_sampling/over_sample_util.py @@ -110,23 +110,16 @@ def slim_index_for_sub_slim_index_via_mask_2d_from( slim_index_for_sub_slim_index = slim_index_for_sub_slim_index_via_mask_2d_from(mask_2d=mask_2d, sub_size=2) """ - total_sub_pixels = total_sub_pixels_2d_from(sub_size=sub_size) + # Step 1: Identify unmasked (False) pixels + unmasked_indices = np.argwhere(~mask_2d) + n_unmasked = unmasked_indices.shape[0] - slim_index_for_sub_slim_index = np.zeros(shape=total_sub_pixels) - slim_index = 0 - sub_slim_index = 0 + # Step 2: Compute total number of sub-pixels + sub_pixels_per_pixel = sub_size ** 2 - for y in range(mask_2d.shape[0]): - for x in range(mask_2d.shape[1]): - if not mask_2d[y, x]: - sub = sub_size[slim_index] - - for y1 in range(sub): - for x1 in range(sub): - slim_index_for_sub_slim_index[sub_slim_index] = slim_index - sub_slim_index += 1 - - slim_index += 1 + # Step 3: Repeat slim indices for each sub-pixel + slim_indices = np.arange(n_unmasked) + slim_index_for_sub_slim_index = np.repeat(slim_indices, sub_pixels_per_pixel) return slim_index_for_sub_slim_index From fbdcf80e68ddd6dd7a795f5f8e64787eafe94159 Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Wed, 18 Jun 2025 21:04:05 +0100 Subject: [PATCH 12/19] removed sub_slim_index_for_sub_native_index_from --- .../over_sampling/over_sample_util.py | 51 ------------ .../over_sample/test_over_sample_util.py | 78 ------------------- 2 files changed, 129 deletions(-) diff --git a/autoarray/operators/over_sampling/over_sample_util.py b/autoarray/operators/over_sampling/over_sample_util.py index b86062d14..24dae44e3 100644 --- a/autoarray/operators/over_sampling/over_sample_util.py +++ b/autoarray/operators/over_sampling/over_sample_util.py @@ -124,57 +124,6 @@ def slim_index_for_sub_slim_index_via_mask_2d_from( return slim_index_for_sub_slim_index -@numba_util.jit() -def sub_slim_index_for_sub_native_index_from(sub_mask_2d: np.ndarray): - """ - Returns a 2D array which maps every `False` entry of a 2D mask to its sub slim mask array. Every - True entry is given a value -1. - - This is used as a convenience tool for creating structures util between different grids and structures. - - For example, if we had a 3x4 mask: - - [[False, True, False, False], - [False, True, False, False], - [False, False, False, True]]] - - The sub_slim_index_for_sub_native_index array would be: - - [[0, -1, 2, 3], - [4, -1, 5, 6], - [7, 8, 9, -1]] - - Parameters - ---------- - sub_mask_2d - The 2D mask that the util array is created for. - - Returns - ------- - ndarray - The 2D array mapping 2D mask entries to their 1D masked array indexes. - - Examples - -------- - mask = np.full(fill_value=False, shape=(9,9)) - sub_two_to_one = mask_to_mask_1d_index_from(mask=mask) - """ - - sub_slim_index_for_sub_native_index = -1 * np.ones(shape=sub_mask_2d.shape) - - sub_mask_1d_index = 0 - - for sub_mask_y in range(sub_mask_2d.shape[0]): - for sub_mask_x in range(sub_mask_2d.shape[1]): - if sub_mask_2d[sub_mask_y, sub_mask_x] == False: - sub_slim_index_for_sub_native_index[sub_mask_y, sub_mask_x] = ( - sub_mask_1d_index - ) - sub_mask_1d_index += 1 - - return sub_slim_index_for_sub_native_index - - @numba_util.jit() def oversample_mask_2d_from(mask: np.ndarray, sub_size: int) -> np.ndarray: """ diff --git a/test_autoarray/operators/over_sample/test_over_sample_util.py b/test_autoarray/operators/over_sample/test_over_sample_util.py index dc6237dc8..be8610d92 100644 --- a/test_autoarray/operators/over_sample/test_over_sample_util.py +++ b/test_autoarray/operators/over_sample/test_over_sample_util.py @@ -79,84 +79,6 @@ def test__slim_index_for_sub_slim_index_via_mask_2d_from(): ).all() -def test__sub_slim_index_for_sub_native_index_from(): - mask = np.full(fill_value=False, shape=(3, 3)) - - sub_mask_1d_index_for_sub_mask_index = ( - util.over_sample.sub_slim_index_for_sub_native_index_from(sub_mask_2d=mask) - ) - - assert ( - sub_mask_1d_index_for_sub_mask_index - == np.array([[0, 1, 2], [3, 4, 5], [6, 7, 8]]) - ).all() - - mask = np.full(fill_value=False, shape=(2, 3)) - - sub_mask_1d_index_for_sub_mask_index = ( - util.over_sample.sub_slim_index_for_sub_native_index_from(sub_mask_2d=mask) - ) - - assert ( - sub_mask_1d_index_for_sub_mask_index == np.array([[0, 1, 2], [3, 4, 5]]) - ).all() - - mask = np.full(fill_value=False, shape=(3, 2)) - - sub_mask_1d_index_for_sub_mask_index = ( - util.over_sample.sub_slim_index_for_sub_native_index_from(sub_mask_2d=mask) - ) - - assert ( - sub_mask_1d_index_for_sub_mask_index == np.array([[0, 1], [2, 3], [4, 5]]) - ).all() - - mask = np.array([[False, True, False], [True, True, False], [False, False, True]]) - - sub_mask_1d_index_for_sub_mask_index = ( - util.over_sample.sub_slim_index_for_sub_native_index_from(sub_mask_2d=mask) - ) - - assert ( - sub_mask_1d_index_for_sub_mask_index - == np.array([[0, -1, 1], [-1, -1, 2], [3, 4, -1]]) - ).all() - - mask = np.array( - [ - [False, True, True, False], - [True, True, False, False], - [False, False, True, False], - ] - ) - - sub_mask_1d_index_for_sub_mask_index = ( - util.over_sample.sub_slim_index_for_sub_native_index_from(sub_mask_2d=mask) - ) - - assert ( - sub_mask_1d_index_for_sub_mask_index - == np.array([[0, -1, -1, 1], [-1, -1, 2, 3], [4, 5, -1, 6]]) - ).all() - - mask = np.array( - [ - [False, True, False], - [True, True, False], - [False, False, True], - [False, False, True], - ] - ) - - sub_mask_1d_index_for_sub_mask_index = ( - util.over_sample.sub_slim_index_for_sub_native_index_from(sub_mask_2d=mask) - ) - - assert ( - sub_mask_1d_index_for_sub_mask_index - == np.array([[0, -1, 1], [-1, -1, 2], [3, 4, -1], [5, 6, -1]]) - ).all() - def test__oversample_mask_from(): mask = np.array( From 91e880502e6fb65a5e5e191b463607cc7b97aaa6 Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Wed, 18 Jun 2025 21:05:56 +0100 Subject: [PATCH 13/19] remove oversample_mask_2d_from --- .../over_sampling/over_sample_util.py | 53 ------------------- .../operators/over_sample/test_decorator.py | 17 +++++- .../over_sample/test_over_sample_util.py | 30 ----------- 3 files changed, 16 insertions(+), 84 deletions(-) diff --git a/autoarray/operators/over_sampling/over_sample_util.py b/autoarray/operators/over_sampling/over_sample_util.py index 24dae44e3..1321ce349 100644 --- a/autoarray/operators/over_sampling/over_sample_util.py +++ b/autoarray/operators/over_sampling/over_sample_util.py @@ -123,59 +123,6 @@ def slim_index_for_sub_slim_index_via_mask_2d_from( return slim_index_for_sub_slim_index - -@numba_util.jit() -def oversample_mask_2d_from(mask: np.ndarray, sub_size: int) -> np.ndarray: - """ - Returns a new mask of shape (mask.shape[0] * sub_size, mask.shape[1] * sub_size) where all boolean values are - expanded according to the `sub_size`. - - For example, if the input mask is: - - mask = np.array([ - [True, True, True], - [True, False, True], - [True, True, True] - ]) - - and the sub_size is 2, the output mask would be: - - expanded_mask = np.array([ - [True, True, True, True, True, True], - [True, True, True, True, True, True], - [True, True, False, False, True, True], - [True, True, False, False, True, True], - [True, True, True, True, True, True], - [True, True, True, True, True, True] - ]) - - This is used throughout the code to handle uniform oversampling calculations. - - Parameters - ---------- - mask - The mask from which the over sample mask is computed. - sub_size - The factor by which the mask is oversampled. - - Returns - ------- - The mask oversampled by the input sub_size. - """ - oversample_mask = np.full( - (mask.shape[0] * sub_size, mask.shape[1] * sub_size), True - ) - - for y in range(mask.shape[0]): - for x in range(mask.shape[1]): - if not mask[y, x]: - oversample_mask[ - y * sub_size : (y + 1) * sub_size, x * sub_size : (x + 1) * sub_size - ] = False - - return oversample_mask - - @numba_util.jit() def sub_size_radial_bins_from( radial_grid: np.ndarray, diff --git a/test_autoarray/operators/over_sample/test_decorator.py b/test_autoarray/operators/over_sample/test_decorator.py index 47dae1e88..42c5054ca 100644 --- a/test_autoarray/operators/over_sample/test_decorator.py +++ b/test_autoarray/operators/over_sample/test_decorator.py @@ -28,7 +28,22 @@ def test__in_grid_2d__over_sample_uniform__out_ndarray_1d(): over_sample_uniform = aa.OverSampler(mask=mask, sub_size=2) - mask_sub_2 = aa.util.over_sample.oversample_mask_2d_from( + def oversample_mask_2d_from(mask: np.ndarray, sub_size: int) -> np.ndarray: + + oversample_mask = np.full( + (mask.shape[0] * sub_size, mask.shape[1] * sub_size), True + ) + + for y in range(mask.shape[0]): + for x in range(mask.shape[1]): + if not mask[y, x]: + oversample_mask[ + y * sub_size: (y + 1) * sub_size, x * sub_size: (x + 1) * sub_size + ] = False + + return oversample_mask + + mask_sub_2 = oversample_mask_2d_from( mask=np.array(mask), sub_size=2 ) diff --git a/test_autoarray/operators/over_sample/test_over_sample_util.py b/test_autoarray/operators/over_sample/test_over_sample_util.py index be8610d92..02bbd6be4 100644 --- a/test_autoarray/operators/over_sample/test_over_sample_util.py +++ b/test_autoarray/operators/over_sample/test_over_sample_util.py @@ -79,36 +79,6 @@ def test__slim_index_for_sub_slim_index_via_mask_2d_from(): ).all() - -def test__oversample_mask_from(): - mask = np.array( - [ - [True, True, True, True], - [True, False, False, True], - [True, False, False, True], - [True, True, True, True], - ] - ) - - oversample_mask = util.over_sample.oversample_mask_2d_from(mask=mask, sub_size=2) - - assert ( - oversample_mask - == np.array( - [ - [True, True, True, True, True, True, True, True], - [True, True, True, True, True, True, True, True], - [True, True, False, False, False, False, True, True], - [True, True, False, False, False, False, True, True], - [True, True, False, False, False, False, True, True], - [True, True, False, False, False, False, True, True], - [True, True, True, True, True, True, True, True], - [True, True, True, True, True, True, True, True], - ] - ) - ).all() - - def test__grid_2d_slim_over_sampled_via_mask_from(): mask = np.array([[True, True, False], [False, False, False], [True, True, False]]) From c4a6339b5f80ce0b9133225661fec3982eeb87fd Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Wed, 18 Jun 2025 21:09:22 +0100 Subject: [PATCH 14/19] sub_size_radial_bins_from no longer uses numba --- .../operators/over_sampling/over_sample_util.py | 17 +++++------------ 1 file changed, 5 insertions(+), 12 deletions(-) diff --git a/autoarray/operators/over_sampling/over_sample_util.py b/autoarray/operators/over_sampling/over_sample_util.py index 1321ce349..a2bd8616d 100644 --- a/autoarray/operators/over_sampling/over_sample_util.py +++ b/autoarray/operators/over_sampling/over_sample_util.py @@ -123,7 +123,6 @@ def slim_index_for_sub_slim_index_via_mask_2d_from( return slim_index_for_sub_slim_index -@numba_util.jit() def sub_size_radial_bins_from( radial_grid: np.ndarray, sub_size_list: np.ndarray, @@ -161,19 +160,13 @@ def sub_size_radial_bins_from( the centre of the mask. """ - sub_size = sub_size_list[-1] * np.ones(radial_grid.shape) + # Use np.searchsorted to find the first index where radial_grid[i] < radial_list[j] + bin_indices = np.searchsorted(radial_list, radial_grid, side="left") - for i in range(radial_grid.shape[0]): - for j in range(len(radial_list)): - if radial_grid[i] < radial_list[j]: - # if use_jax: - # # while this makes it run, it is very, very slow - # sub_size = sub_size.at[i].set(sub_size_list[j]) - # else: - sub_size[i] = sub_size_list[j] - break + # Clip indices to stay within bounds of sub_size_list + bin_indices = np.clip(bin_indices, 0, len(sub_size_list) - 1) - return sub_size + return sub_size_list[bin_indices] @numba_util.jit() From 7b1c3b2de6540827efb9ed9583da32809d2eb53e Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Thu, 19 Jun 2025 14:34:39 +0100 Subject: [PATCH 15/19] grid_2d_slim_over_sampled_via_mask_from does not use mask --- .../over_sampling/over_sample_util.py | 169 ++++++------------ .../operators/over_sampling/over_sampler.py | 14 ++ .../over_sample/test_over_sample_util.py | 2 +- 3 files changed, 68 insertions(+), 117 deletions(-) diff --git a/autoarray/operators/over_sampling/over_sample_util.py b/autoarray/operators/over_sampling/over_sample_util.py index a2bd8616d..fba8e8a51 100644 --- a/autoarray/operators/over_sampling/over_sample_util.py +++ b/autoarray/operators/over_sampling/over_sample_util.py @@ -1,4 +1,5 @@ from __future__ import annotations +from collections import defaultdict import numpy as np from typing import TYPE_CHECKING, Union from typing import List, Tuple @@ -169,7 +170,6 @@ def sub_size_radial_bins_from( return sub_size_list[bin_indices] -@numba_util.jit() def grid_2d_slim_over_sampled_via_mask_from( mask_2d: np.ndarray, pixel_scales: ty.PixelScales, @@ -215,121 +215,58 @@ def grid_2d_slim_over_sampled_via_mask_from( grid_slim = grid_2d_slim_over_sampled_via_mask_from(mask=mask, pixel_scales=(0.5, 0.5), sub_size=1, origin=(0.0, 0.0)) """ - total_sub_pixels = np.sum(sub_size**2) - - grid_slim = np.zeros(shape=(total_sub_pixels, 2)) - - centres_scaled = geometry_util.central_scaled_coordinate_2d_numba_from( - shape_native=mask_2d.shape, pixel_scales=pixel_scales, origin=origin - ) - - index = 0 - sub_index = 0 - - for y in range(mask_2d.shape[0]): - for x in range(mask_2d.shape[1]): - if not mask_2d[y, x]: - sub = sub_size[index] - - y_sub_half = pixel_scales[0] / 2 - y_sub_step = pixel_scales[0] / (sub) - - x_sub_half = pixel_scales[1] / 2 - x_sub_step = pixel_scales[1] / (sub) - - y_scaled = (y - centres_scaled[0]) * pixel_scales[0] - x_scaled = (x - centres_scaled[1]) * pixel_scales[1] - - for y1 in range(sub): - for x1 in range(sub): - grid_slim[sub_index, 0] = -( - y_scaled - y_sub_half + y1 * y_sub_step + (y_sub_step / 2.0) - ) - grid_slim[sub_index, 1] = ( - x_scaled - x_sub_half + x1 * x_sub_step + (x_sub_step / 2.0) - ) - sub_index += 1 - - index += 1 - - return grid_slim - - -@numba_util.jit() -def binned_array_2d_from( - array_2d: np.ndarray, - mask_2d: np.ndarray, - sub_size: np.ndarray, -) -> np.ndarray: - """ - For a sub-grid, every unmasked pixel of its 2D mask with shape (total_y_pixels, total_x_pixels) is divided into - a finer uniform grid of shape (total_y_pixels*sub_size, total_x_pixels*sub_size). This routine computes the (y,x) - scaled coordinates a the centre of every sub-pixel defined by this 2D mask array. - - The sub-grid is returned on an array of shape (total_unmasked_pixels*sub_size**2, 2). y coordinates are - stored in the 0 index of the second dimension, x coordinates in the 1 index. Masked coordinates are therefore - removed and not included in the slimmed grid. - - Grid2D are defined from the top-left corner, where the first unmasked sub-pixel corresponds to index 0. - Sub-pixels that are part of the same mask array pixel are indexed next to one another, such that the second - sub-pixel in the first pixel has index 1, its next sub-pixel has index 2, and so forth. - - Parameters - ---------- - mask_2d - A 2D array of bools, where `False` values are unmasked and therefore included as part of the calculated - sub-grid. - pixel_scales - The (y,x) scaled units to pixel units conversion factor of the 2D mask array. - sub_size - The size of the sub-grid that each pixel of the 2D mask array is divided into. - origin - The (y,x) origin of the 2D array, which the sub-grid is shifted around. - - Returns - ------- - ndarray - A slimmed sub grid of (y,x) scaled coordinates at the centre of every pixel unmasked pixel on the 2D mask - array. The sub grid array has dimensions (total_unmasked_pixels*sub_size**2, 2). - - Examples - -------- - mask = np.array([[True, False, True], - [False, False, False] - [True, False, True]]) - grid_slim = grid_2d_slim_over_sampled_via_mask_from(mask=mask, pixel_scales=(0.5, 0.5), sub_size=1, origin=(0.0, 0.0)) - """ - - total_pixels = np.sum(~mask_2d) - - sub_fraction = 1.0 / sub_size**2 - - binned_array_2d_slim = np.zeros(shape=total_pixels) - - index = 0 - sub_index = 0 - - for y in range(mask_2d.shape[0]): - for x in range(mask_2d.shape[1]): - if not mask_2d[y, x]: - sub = sub_size[index] - - for y1 in range(sub): - for x1 in range(sub): - # if use_jax: - # binned_array_2d_slim = binned_array_2d_slim.at[index].add( - # array_2d[sub_index] * sub_fraction[index] - # ) - # else: - binned_array_2d_slim[index] += ( - array_2d[sub_index] * sub_fraction[index] - ) - sub_index += 1 - - index += 1 - - return binned_array_2d_slim - + H, W = mask_2d.shape + sy, sx = pixel_scales + oy, ox = origin + + # 1) Find unmasked pixels in row-major order + rows, cols = np.nonzero(~mask_2d) + Npix = rows.size + + # 2) Normalize sub_size input + sub_arr = np.asarray(sub_size) + sub_arr = np.full(Npix, sub_arr, dtype=int) if sub_arr.size == 1 else sub_arr + + # 3) Pixel centers in physical coords, y↑up + cy = (H - 1) / 2.0 + cx = (W - 1) / 2.0 + y_pix = (cy - rows) * sy + oy + x_pix = (cols - cx) * sx + ox + + # Pre‐group pixel indices by sub_size + groups = defaultdict(list) + for i, s in enumerate(sub_arr): + groups[s].append(i) + + # Prepare output + total = np.sum(sub_arr * sub_arr) + coords = np.empty((total, 2), float) + idx = 0 + + for s, pix_indices in groups.items(): + # Compute offsets once for this sub_size + dy, dx = sy/s, sx/s + y_off = np.linspace(+sy/2 - dy/2, -sy/2 + dy/2, s) + x_off = np.linspace(-sx/2 + dx/2, +sx/2 - dx/2, s) + y_sub, x_sub = np.meshgrid(y_off, x_off, indexing="ij") + y_sub = y_sub.ravel() + x_sub = x_sub.ravel() + n_sub = s*s + + # Now vectorize over all pixels in this group + pix_idx = np.array(pix_indices) + y_centers = y_pix[pix_idx] + x_centers = x_pix[pix_idx] + + # Repeat‐tile to shape (len(pix_idx)*n_sub,) + all_y = np.repeat(y_centers, n_sub) + np.tile(y_sub, len(pix_idx)) + all_x = np.repeat(x_centers, n_sub) + np.tile(x_sub, len(pix_idx)) + + coords[idx:idx + all_y.size, 0] = all_y + coords[idx:idx + all_x.size, 1] = all_x + idx += all_y.size + + return coords def over_sample_size_via_radial_bins_from( grid: Grid2D, diff --git a/autoarray/operators/over_sampling/over_sampler.py b/autoarray/operators/over_sampling/over_sampler.py index 268cc7be9..a71389eff 100644 --- a/autoarray/operators/over_sampling/over_sampler.py +++ b/autoarray/operators/over_sampling/over_sampler.py @@ -230,6 +230,20 @@ def binned_array_2d_from(self, array: Array2D) -> "Array2D": In **PyAutoCTI** all `Array2D` objects are used in their `native` representation without sub-gridding. Significant memory can be saved by only store this format, thus the `native_binned_only` config override can force this behaviour. It is recommended users do not use this option to avoid unexpected behaviour. + + Old docstring: + + For a sub-grid, every unmasked pixel of its 2D mask with shape (total_y_pixels, total_x_pixels) is divided into + a finer uniform grid of shape (total_y_pixels*sub_size, total_x_pixels*sub_size). This routine computes the (y,x) + scaled coordinates a the centre of every sub-pixel defined by this 2D mask array. + + The sub-grid is returned on an array of shape (total_unmasked_pixels*sub_size**2, 2). y coordinates are + stored in the 0 index of the second dimension, x coordinates in the 1 index. Masked coordinates are therefore + removed and not included in the slimmed grid. + + Grid2D are defined from the top-left corner, where the first unmasked sub-pixel corresponds to index 0. + Sub-pixels that are part of the same mask array pixel are indexed next to one another, such that the second + sub-pixel in the first pixel has index 1, its next sub-pixel has index 2, and so forth. """ if conf.instance["general"]["structures"]["native_binned_only"]: return self diff --git a/test_autoarray/operators/over_sample/test_over_sample_util.py b/test_autoarray/operators/over_sample/test_over_sample_util.py index 02bbd6be4..fb68da898 100644 --- a/test_autoarray/operators/over_sample/test_over_sample_util.py +++ b/test_autoarray/operators/over_sample/test_over_sample_util.py @@ -83,7 +83,7 @@ def test__grid_2d_slim_over_sampled_via_mask_from(): mask = np.array([[True, True, False], [False, False, False], [True, True, False]]) grid = aa.util.over_sample.grid_2d_slim_over_sampled_via_mask_from( - mask_2d=mask, pixel_scales=(3.0, 3.0), sub_size=np.array([2, 2, 2, 2, 2]) + mask_2d=mask, pixel_scales=(3.0, 3.0), sub_size=2 ) assert ( From d5d28c968a07d9005bd669f6a1ae67871c321a1f Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Thu, 19 Jun 2025 14:39:55 +0100 Subject: [PATCH 16/19] removed final numba functions not used for inversion --- autoarray/geometry/geometry_util.py | 71 +------------------ .../over_sampling/over_sample_util.py | 18 ++--- .../operators/over_sampling/over_sampler.py | 1 - autoarray/operators/transformer.py | 26 +++---- autoarray/structures/arrays/array_2d_util.py | 10 ++- autoarray/structures/grids/grid_2d_util.py | 1 + .../operators/over_sample/test_decorator.py | 7 +- .../structures/arrays/test_array_2d_util.py | 4 +- .../structures/grids/test_grid_2d_util.py | 2 +- 9 files changed, 38 insertions(+), 102 deletions(-) diff --git a/autoarray/geometry/geometry_util.py b/autoarray/geometry/geometry_util.py index 4cf32a082..16ee5c460 100644 --- a/autoarray/geometry/geometry_util.py +++ b/autoarray/geometry/geometry_util.py @@ -2,8 +2,6 @@ import numpy as np from typing import Tuple, Union - -from autoarray import numba_util from autoarray import type as ty @@ -178,70 +176,6 @@ def convert_pixel_scales_2d(pixel_scales: ty.PixelScales) -> Tuple[float, float] return pixel_scales -@numba_util.jit() -def central_pixel_coordinates_2d_numba_from( - shape_native: Tuple[int, int], -) -> Tuple[float, float]: - """ - Returns the central pixel coordinates of a 2D geometry (and therefore a 2D data structure like an ``Array2D``) - from the shape of that data structure. - - Examples of the central pixels are as follows: - - - For a 3x3 image, the central pixel is pixel [1, 1]. - - For a 4x4 image, the central pixel is [1.5, 1.5]. - - Parameters - ---------- - shape_native - The dimensions of the data structure, which can be in 1D, 2D or higher dimensions. - - Returns - ------- - The central pixel coordinates of the data structure. - """ - return (float(shape_native[0] - 1) / 2, float(shape_native[1] - 1) / 2) - - -@numba_util.jit() -def central_scaled_coordinate_2d_numba_from( - shape_native: Tuple[int, int], - pixel_scales: ty.PixelScales, - origin: Tuple[float, float] = (0.0, 0.0), -) -> Tuple[float, float]: - """ - Returns the central scaled coordinates of a 2D geometry (and therefore a 2D data structure like an ``Array2D``) - from the shape of that data structure. - - This is computed by using the data structure's shape and converting it to scaled units using an input - pixel-coordinates to scaled-coordinate conversion factor `pixel_scales`. - - The origin of the scaled grid can also be input and moved from (0.0, 0.0). - - Parameters - ---------- - shape_native - The 2D shape of the data structure whose central scaled coordinates are computed. - pixel_scales - The (y,x) scaled units to pixel units conversion factor of the 2D data structure. - origin - The (y,x) scaled units origin of the coordinate system the central scaled coordinate is computed on. - - Returns - ------- - The central coordinates of the 2D data structure in scaled units. - """ - - central_pixel_coordinates = central_pixel_coordinates_2d_numba_from( - shape_native=shape_native - ) - - y_pixel = central_pixel_coordinates[0] + (origin[0] / pixel_scales[0]) - x_pixel = central_pixel_coordinates[1] - (origin[1] / pixel_scales[1]) - - return (y_pixel, x_pixel) - - def central_pixel_coordinates_2d_from( shape_native: Tuple[int, int], ) -> Tuple[float, float]: @@ -294,7 +228,7 @@ def central_scaled_coordinate_2d_from( The central coordinates of the 2D data structure in scaled units. """ - central_pixel_coordinates = central_pixel_coordinates_2d_numba_from( + central_pixel_coordinates = central_pixel_coordinates_2d_from( shape_native=shape_native ) @@ -367,7 +301,6 @@ def pixel_coordinates_2d_from( return (y_pixel, x_pixel) -@numba_util.jit() def scaled_coordinates_2d_from( pixel_coordinates_2d: Tuple[float, float], shape_native: Tuple[int, int], @@ -411,7 +344,7 @@ def scaled_coordinates_2d_from( origin=(0.0, 0.0) ) """ - central_scaled_coordinates = central_scaled_coordinate_2d_numba_from( + central_scaled_coordinates = central_scaled_coordinate_2d_from( shape_native=shape_native, pixel_scales=pixel_scales, origin=origins ) diff --git a/autoarray/operators/over_sampling/over_sample_util.py b/autoarray/operators/over_sampling/over_sample_util.py index fba8e8a51..c038a5491 100644 --- a/autoarray/operators/over_sampling/over_sample_util.py +++ b/autoarray/operators/over_sampling/over_sample_util.py @@ -9,7 +9,6 @@ if TYPE_CHECKING: from autoarray.structures.grids.uniform_2d import Grid2D -from autoarray.geometry import geometry_util from autoarray.mask.mask_2d import Mask2D from autoarray import numba_util @@ -50,7 +49,6 @@ def over_sample_size_convert_to_array_2d_from( return Array2D(values=np.array(over_sample_size).astype("int"), mask=mask) -@numba_util.jit() def total_sub_pixels_2d_from(sub_size: np.ndarray) -> int: """ Returns the total number of sub-pixels in unmasked pixels in a mask. @@ -116,7 +114,7 @@ def slim_index_for_sub_slim_index_via_mask_2d_from( n_unmasked = unmasked_indices.shape[0] # Step 2: Compute total number of sub-pixels - sub_pixels_per_pixel = sub_size ** 2 + sub_pixels_per_pixel = sub_size**2 # Step 3: Repeat slim indices for each sub-pixel slim_indices = np.arange(n_unmasked) @@ -124,6 +122,7 @@ def slim_index_for_sub_slim_index_via_mask_2d_from( return slim_index_for_sub_slim_index + def sub_size_radial_bins_from( radial_grid: np.ndarray, sub_size_list: np.ndarray, @@ -245,13 +244,13 @@ def grid_2d_slim_over_sampled_via_mask_from( for s, pix_indices in groups.items(): # Compute offsets once for this sub_size - dy, dx = sy/s, sx/s - y_off = np.linspace(+sy/2 - dy/2, -sy/2 + dy/2, s) - x_off = np.linspace(-sx/2 + dx/2, +sx/2 - dx/2, s) + dy, dx = sy / s, sx / s + y_off = np.linspace(+sy / 2 - dy / 2, -sy / 2 + dy / 2, s) + x_off = np.linspace(-sx / 2 + dx / 2, +sx / 2 - dx / 2, s) y_sub, x_sub = np.meshgrid(y_off, x_off, indexing="ij") y_sub = y_sub.ravel() x_sub = x_sub.ravel() - n_sub = s*s + n_sub = s * s # Now vectorize over all pixels in this group pix_idx = np.array(pix_indices) @@ -262,12 +261,13 @@ def grid_2d_slim_over_sampled_via_mask_from( all_y = np.repeat(y_centers, n_sub) + np.tile(y_sub, len(pix_idx)) all_x = np.repeat(x_centers, n_sub) + np.tile(x_sub, len(pix_idx)) - coords[idx:idx + all_y.size, 0] = all_y - coords[idx:idx + all_x.size, 1] = all_x + coords[idx : idx + all_y.size, 0] = all_y + coords[idx : idx + all_x.size, 1] = all_x idx += all_y.size return coords + def over_sample_size_via_radial_bins_from( grid: Grid2D, sub_size_list: List[int], diff --git a/autoarray/operators/over_sampling/over_sampler.py b/autoarray/operators/over_sampling/over_sampler.py index a71389eff..29c93aa7e 100644 --- a/autoarray/operators/over_sampling/over_sampler.py +++ b/autoarray/operators/over_sampling/over_sampler.py @@ -272,7 +272,6 @@ def binned_array_2d_from(self, array: Array2D) -> "Array2D": mask=self.mask, ) - @cached_property def slim_for_sub_slim(self) -> np.ndarray: """ diff --git a/autoarray/operators/transformer.py b/autoarray/operators/transformer.py index 675123e78..ad861f0fb 100644 --- a/autoarray/operators/transformer.py +++ b/autoarray/operators/transformer.py @@ -418,22 +418,22 @@ def image_from( def transform_mapping_matrix(self, mapping_matrix: np.ndarray) -> np.ndarray: """ - Applies the NUFFT forward transform to each column of a mapping matrix, producing transformed visibilities. + Applies the NUFFT forward transform to each column of a mapping matrix, producing transformed visibilities. - Parameters - ---------- - mapping_matrix - A 2D array where each column corresponds to a source-plane pixel intensity distribution flattened into image space. + Parameters + ---------- + mapping_matrix + A 2D array where each column corresponds to a source-plane pixel intensity distribution flattened into image space. - Returns - ------- - A complex-valued 2D array where each column contains the visibilities corresponding to the respective column - in the input mapping matrix. + Returns + ------- + A complex-valued 2D array where each column contains the visibilities corresponding to the respective column + in the input mapping matrix. - Notes - ----- - - Each column of the input mapping matrix is reshaped into the native 2D image grid before transformation. - - This method repeatedly calls `visibilities_from` for each column, which may be computationally intensive. + Notes + ----- + - Each column of the input mapping matrix is reshaped into the native 2D image grid before transformation. + - This method repeatedly calls `visibilities_from` for each column, which may be computationally intensive. """ transformed_mapping_matrix = 0 + 0j * np.zeros( (self.uv_wavelengths.shape[0], mapping_matrix.shape[1]) diff --git a/autoarray/structures/arrays/array_2d_util.py b/autoarray/structures/arrays/array_2d_util.py index 92e99d43a..aff0872b5 100644 --- a/autoarray/structures/arrays/array_2d_util.py +++ b/autoarray/structures/arrays/array_2d_util.py @@ -337,10 +337,13 @@ def resized_array_2d_from( dst_x_end = dst_x_start + (src_x_end - src_x_start) # Copy overlapping region from source to destination - resized_array[dst_y_start:dst_y_end, dst_x_start:dst_x_end] = array_2d[src_y_start:src_y_end, src_x_start:src_x_end] + resized_array[dst_y_start:dst_y_end, dst_x_start:dst_x_end] = array_2d[ + src_y_start:src_y_end, src_x_start:src_x_end + ] return resized_array + def index_2d_for_index_slim_from(indexes_slim: np.ndarray, shape_native) -> np.ndarray: """ For pixels on a native 2D array of shape (total_y_pixels, total_x_pixels), this array maps the slimmed 1D pixel @@ -384,6 +387,7 @@ def index_2d_for_index_slim_from(indexes_slim: np.ndarray, shape_native) -> np.n return index_2d_for_index_slim + def index_slim_for_index_2d_from(indexes_2d: np.ndarray, shape_native) -> np.ndarray: """ For pixels on a native 2D array of shape (total_y_pixels, total_x_pixels), this array maps the 2D pixel indexes to @@ -417,7 +421,9 @@ def index_slim_for_index_2d_from(indexes_2d: np.ndarray, shape_native) -> np.nda indexes_flat = index_flat_for_index_2d_from(indexes_2d=indexes_2d, shape=(3,3)) """ # Calculate 1D indexes as row_index * number_of_columns + col_index - index_slim_for_index_native_2d = indexes_2d[:, 0] * shape_native[1] + indexes_2d[:, 1] + index_slim_for_index_native_2d = ( + indexes_2d[:, 0] * shape_native[1] + indexes_2d[:, 1] + ) return index_slim_for_index_native_2d diff --git a/autoarray/structures/grids/grid_2d_util.py b/autoarray/structures/grids/grid_2d_util.py index 00c67cd7d..ea298bce2 100644 --- a/autoarray/structures/grids/grid_2d_util.py +++ b/autoarray/structures/grids/grid_2d_util.py @@ -395,6 +395,7 @@ def grid_2d_via_shape_native_from( origin=origin, ) + def _radial_projected_shape_slim_from( extent: np.ndarray, centre: Tuple[float, float], diff --git a/test_autoarray/operators/over_sample/test_decorator.py b/test_autoarray/operators/over_sample/test_decorator.py index 42c5054ca..72af0215a 100644 --- a/test_autoarray/operators/over_sample/test_decorator.py +++ b/test_autoarray/operators/over_sample/test_decorator.py @@ -38,14 +38,13 @@ def oversample_mask_2d_from(mask: np.ndarray, sub_size: int) -> np.ndarray: for x in range(mask.shape[1]): if not mask[y, x]: oversample_mask[ - y * sub_size: (y + 1) * sub_size, x * sub_size: (x + 1) * sub_size + y * sub_size : (y + 1) * sub_size, + x * sub_size : (x + 1) * sub_size, ] = False return oversample_mask - mask_sub_2 = oversample_mask_2d_from( - mask=np.array(mask), sub_size=2 - ) + mask_sub_2 = oversample_mask_2d_from(mask=np.array(mask), sub_size=2) mask_sub_2 = aa.Mask2D(mask=mask_sub_2, pixel_scales=(0.5, 0.5)) diff --git a/test_autoarray/structures/arrays/test_array_2d_util.py b/test_autoarray/structures/arrays/test_array_2d_util.py index 8b8056f3b..0c6ae8bc7 100644 --- a/test_autoarray/structures/arrays/test_array_2d_util.py +++ b/test_autoarray/structures/arrays/test_array_2d_util.py @@ -247,6 +247,7 @@ def test__resized_array_2d_from__padding(): ) ).all() + def test__resized_array_2d_from__padding_with_new_origin(): array = np.ones((3, 3)) array[2, 2] = 2.0 @@ -286,7 +287,6 @@ def test__resized_array_2d_from__padding_with_new_origin(): ).all() - def test__index_2d_for_index_slim_from(): indexes_1d = np.array([0, 1, 2, 3, 4, 5, 6, 7, 8]) @@ -431,7 +431,6 @@ def test__array_2d_slim_from(): assert (array_2d_slim == np.array([2, 4, 5, 6, 8])).all() - def test__array_2d_native_from(): array_2d_slim = np.array([1.0, 2.0, 3.0, 4.0]) @@ -476,4 +475,3 @@ def test__array_2d_native_from(): [[1.0, 2.0, 0.0, 0.0], [3.0, 0.0, 0.0, 0.0], [-1.0, -2.0, 0.0, -3.0]] ) ).all() - diff --git a/test_autoarray/structures/grids/test_grid_2d_util.py b/test_autoarray/structures/grids/test_grid_2d_util.py index 4adc6c020..034f267e5 100644 --- a/test_autoarray/structures/grids/test_grid_2d_util.py +++ b/test_autoarray/structures/grids/test_grid_2d_util.py @@ -460,4 +460,4 @@ def test__grid_2d_native_from(): [[-1.0, -1.0], [-2.0, -2.0], [0.0, 0.0], [-3.0, -3.0]], ] ) - ).all() \ No newline at end of file + ).all() From 8fca9ba9b5e208490261a56dd067a8b430ad3154 Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Thu, 19 Jun 2025 14:40:45 +0100 Subject: [PATCH 17/19] numba stops pixelization if not installed --- autoarray/inversion/inversion/abstract.py | 12 ++++++++++++ 1 file changed, 12 insertions(+) diff --git a/autoarray/inversion/inversion/abstract.py b/autoarray/inversion/inversion/abstract.py index bec0b1ce2..daecbd4e2 100644 --- a/autoarray/inversion/inversion/abstract.py +++ b/autoarray/inversion/inversion/abstract.py @@ -74,6 +74,18 @@ def __init__( A dictionary which contains timing of certain functions calls which is used for profiling. """ + try: + import numba + except ModuleNotFoundError: + raise exc.InversionException( + "Inversion functionality (linear light profiles, pixelized reconstructions) is " + "disabled if numba is not installed.\n\n" + "This is because the run-times without numba are too slow.\n\n" + "Please install numba, which is described at the following web page:\n\n" + "https://pyautolens.readthedocs.io/en/latest/installation/overview.html" + ) + + self.dataset = dataset self.linear_obj_list = linear_obj_list From 5c41ba11db725595302e3cf085b0b7e5cf17ddf1 Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Thu, 19 Jun 2025 16:18:16 +0100 Subject: [PATCH 18/19] test --- autoarray/operators/over_sampling/over_sample_util.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/autoarray/operators/over_sampling/over_sample_util.py b/autoarray/operators/over_sampling/over_sample_util.py index a2bd8616d..653c55ed3 100644 --- a/autoarray/operators/over_sampling/over_sample_util.py +++ b/autoarray/operators/over_sampling/over_sample_util.py @@ -169,7 +169,7 @@ def sub_size_radial_bins_from( return sub_size_list[bin_indices] -@numba_util.jit() +# @numba_util.jit() def grid_2d_slim_over_sampled_via_mask_from( mask_2d: np.ndarray, pixel_scales: ty.PixelScales, From fdac408d5991db18a8d146dcdcd5e3f7aa54ffce Mon Sep 17 00:00:00 2001 From: Jammy2211 Date: Thu, 19 Jun 2025 16:40:15 +0100 Subject: [PATCH 19/19] use simpler grid_2d_slim_over_sampled_via_mask_from which works --- .../over_sampling/over_sample_util.py | 66 +++++++++---------- 1 file changed, 30 insertions(+), 36 deletions(-) diff --git a/autoarray/operators/over_sampling/over_sample_util.py b/autoarray/operators/over_sampling/over_sample_util.py index c038a5491..e32c377ad 100644 --- a/autoarray/operators/over_sampling/over_sample_util.py +++ b/autoarray/operators/over_sampling/over_sample_util.py @@ -1,5 +1,4 @@ from __future__ import annotations -from collections import defaultdict import numpy as np from typing import TYPE_CHECKING, Union from typing import List, Tuple @@ -11,7 +10,6 @@ from autoarray.mask.mask_2d import Mask2D -from autoarray import numba_util from autoarray import type as ty @@ -168,7 +166,6 @@ def sub_size_radial_bins_from( return sub_size_list[bin_indices] - def grid_2d_slim_over_sampled_via_mask_from( mask_2d: np.ndarray, pixel_scales: ty.PixelScales, @@ -218,54 +215,51 @@ def grid_2d_slim_over_sampled_via_mask_from( sy, sx = pixel_scales oy, ox = origin - # 1) Find unmasked pixels in row-major order + # 1) Find unmasked pixel indices in row-major order rows, cols = np.nonzero(~mask_2d) Npix = rows.size - # 2) Normalize sub_size input + # 2) Broadcast or validate sub_size array sub_arr = np.asarray(sub_size) - sub_arr = np.full(Npix, sub_arr, dtype=int) if sub_arr.size == 1 else sub_arr - - # 3) Pixel centers in physical coords, y↑up + if sub_arr.ndim == 0: + sub_arr = np.full(Npix, int(sub_arr), int) + elif sub_arr.ndim == 1 and sub_arr.size == Npix: + sub_arr = sub_arr.astype(int) + else: + raise ValueError(f"sub_size must be scalar or length-{Npix} array, got shape {sub_arr.shape}") + + # 3) Compute pixel centers (y ↑ up, x → right) cy = (H - 1) / 2.0 cx = (W - 1) / 2.0 y_pix = (cy - rows) * sy + oy x_pix = (cols - cx) * sx + ox - # Pre‐group pixel indices by sub_size - groups = defaultdict(list) - for i, s in enumerate(sub_arr): - groups[s].append(i) - - # Prepare output - total = np.sum(sub_arr * sub_arr) - coords = np.empty((total, 2), float) - idx = 0 - - for s, pix_indices in groups.items(): - # Compute offsets once for this sub_size - dy, dx = sy / s, sx / s - y_off = np.linspace(+sy / 2 - dy / 2, -sy / 2 + dy / 2, s) - x_off = np.linspace(-sx / 2 + dx / 2, +sx / 2 - dx / 2, s) + # 4) For each pixel, generate its sub-pixel coords and collect + coords_list = [] + for i in range(Npix): + s = sub_arr[i] + dy = sy / s + dx = sx / s + + # y offsets: from top (+sy/2 - dy/2) down to bottom (-sy/2 + dy/2) + y_off = np.linspace(+sy/2 - dy/2, -sy/2 + dy/2, s) + # x offsets: left to right + x_off = np.linspace(-sx/2 + dx/2, +sx/2 - dx/2, s) + + # build subgrid y_sub, x_sub = np.meshgrid(y_off, x_off, indexing="ij") y_sub = y_sub.ravel() x_sub = x_sub.ravel() - n_sub = s * s - - # Now vectorize over all pixels in this group - pix_idx = np.array(pix_indices) - y_centers = y_pix[pix_idx] - x_centers = x_pix[pix_idx] - # Repeat‐tile to shape (len(pix_idx)*n_sub,) - all_y = np.repeat(y_centers, n_sub) + np.tile(y_sub, len(pix_idx)) - all_x = np.repeat(x_centers, n_sub) + np.tile(x_sub, len(pix_idx)) + # center + offsets + y_center = y_pix[i] + x_center = x_pix[i] + coords = np.stack([y_center + y_sub, x_center + x_sub], axis=1) - coords[idx : idx + all_y.size, 0] = all_y - coords[idx : idx + all_x.size, 1] = all_x - idx += all_y.size + coords_list.append(coords) - return coords + # 5) Concatenate all sub-pixel blocks in row-major pixel order + return np.vstack(coords_list) def over_sample_size_via_radial_bins_from(