diff --git a/autoarray/dataset/grids.py b/autoarray/dataset/grids.py
index d97fd3f4d..18115d20f 100644
--- a/autoarray/dataset/grids.py
+++ b/autoarray/dataset/grids.py
@@ -17,6 +17,7 @@ def __init__(
         over_sample_size_lp: Union[int, Array2D],
         over_sample_size_pixelization: Union[int, Array2D],
         psf: Optional[Kernel2D] = None,
+        use_w_tilde: bool = False,
     ):
         """
         Contains grids of (y,x) Cartesian coordinates at the centre of every pixel in the dataset's image and
@@ -64,13 +65,11 @@ def __init__(
             mask=self.mask,
             over_sample_size=self.over_sample_size_lp,
         )
-        self.lp.over_sampled
 
         self.pixelization = Grid2D.from_mask(
             mask=self.mask,
             over_sample_size=self.over_sample_size_pixelization,
         )
-        self.pixelization.over_sampled
 
         if self.psf is None:
             self.blurring = None
@@ -79,12 +78,13 @@ def __init__(
                 self.blurring = self.lp.blurring_grid_via_kernel_shape_from(
                     kernel_shape_native=self.psf.shape_native,
                 )
-                self.blurring.over_sampled
             except exc.MaskException:
                 self.blurring = None
 
         self.border_relocator = BorderRelocator(
-            mask=self.mask, sub_size=self.over_sample_size_pixelization
+            mask=self.mask,
+            sub_size=self.over_sample_size_pixelization,
+            use_w_tilde=use_w_tilde,
         )
 
 
diff --git a/autoarray/dataset/imaging/dataset.py b/autoarray/dataset/imaging/dataset.py
index 1eef061fa..a98a443ac 100644
--- a/autoarray/dataset/imaging/dataset.py
+++ b/autoarray/dataset/imaging/dataset.py
@@ -191,11 +191,14 @@ def __init__(
             if psf.mask.shape[0] % 2 == 0 or psf.mask.shape[1] % 2 == 0:
                 raise exc.KernelException("Kernel2D Kernel2D must be odd")
 
+        use_w_tilde = True if w_tilde is not None else False
+
         self.grids = GridsDataset(
             mask=self.data.mask,
             over_sample_size_lp=self.over_sample_size_lp,
             over_sample_size_pixelization=self.over_sample_size_pixelization,
             psf=self.psf,
+            use_w_tilde=use_w_tilde,
         )
 
         self.w_tilde = w_tilde
diff --git a/autoarray/dataset/interferometer/dataset.py b/autoarray/dataset/interferometer/dataset.py
index aa62fb3a2..9af9de286 100644
--- a/autoarray/dataset/interferometer/dataset.py
+++ b/autoarray/dataset/interferometer/dataset.py
@@ -100,10 +100,13 @@ def __init__(
 
         self.dft_preload_transform = dft_preload_transform
 
+        use_w_tilde = True if w_tilde is not None else False
+
         self.grids = GridsDataset(
             mask=self.real_space_mask,
             over_sample_size_lp=self.over_sample_size_lp,
             over_sample_size_pixelization=self.over_sample_size_pixelization,
+            use_w_tilde=use_w_tilde,
         )
 
         self.w_tilde = w_tilde
diff --git a/autoarray/inversion/inversion/abstract.py b/autoarray/inversion/inversion/abstract.py
index 4266cb086..7b809b1d5 100644
--- a/autoarray/inversion/inversion/abstract.py
+++ b/autoarray/inversion/inversion/abstract.py
@@ -293,7 +293,7 @@ def mapping_matrix(self) -> np.ndarray:
     def operated_mapping_matrix_list(self) -> np.ndarray:
         raise NotImplementedError
 
-    @property
+    @cached_property
     def operated_mapping_matrix(self) -> np.ndarray:
         """
         The `operated_mapping_matrix` of a linear object describes the mappings between the observed data's values and
@@ -314,7 +314,7 @@ def data_vector(self) -> np.ndarray:
     def curvature_matrix(self) -> np.ndarray:
         raise NotImplementedError
 
-    @property
+    @cached_property
     def regularization_matrix(self) -> Optional[np.ndarray]:
         """
         The regularization matrix H is used to impose smoothness on our inversion's reconstruction. This enters the
@@ -346,7 +346,7 @@ def regularization_matrix(self) -> Optional[np.ndarray]:
             *[linear_obj.regularization_matrix for linear_obj in self.linear_obj_list]
         )
 
-    @property
+    @cached_property
     def regularization_matrix_reduced(self) -> Optional[np.ndarray]:
         """
         The regularization matrix H is used to impose smoothness on our inversion's reconstruction. This enters the
@@ -360,7 +360,6 @@ def regularization_matrix_reduced(self) -> Optional[np.ndarray]:
         The scipy function `block_diag` has an overhead associated with it and if there is only one mapper and
         regularization it is bypassed.
         """
-
         if self.all_linear_obj_have_regularization:
             return self.regularization_matrix
 
@@ -380,7 +379,6 @@ def curvature_reg_matrix(self) -> np.ndarray:
         to ensure if we access it after computing the `curvature_reg_matrix` it is correctly recalculated in a new
         array of memory.
         """
-
         if not self.has(cls=AbstractRegularization):
             return self.curvature_matrix
 
@@ -400,7 +398,6 @@ def curvature_reg_matrix_reduced(self) -> Optional[np.ndarray]:
         The scipy function `block_diag` has an overhead associated with it and if there is only one mapper and
         regularization it is bypassed.
         """
-
         if self.all_linear_obj_have_regularization:
             return self.curvature_reg_matrix
 
@@ -426,7 +423,6 @@ def reconstruction(self) -> np.ndarray:
             ZTZ := np.dot(Z.T, Z)
             ZTx := np.dot(Z.T, x)
         """
-
         if self.settings.use_positive_only_solver:
 
             if (
@@ -481,7 +477,7 @@ def reconstruction(self) -> np.ndarray:
             xp=self._xp,
         )
 
-    @property
+    @cached_property
     def reconstruction_reduced(self) -> np.ndarray:
         """
         Solve the linear system [F + reg_coeff*H] S = D -> S = [F + reg_coeff*H]^-1 D given by equation (12)
@@ -489,7 +485,6 @@ def reconstruction_reduced(self) -> np.ndarray:
 
         S is the vector of reconstructed inversion values.
         """
-
         if self.all_linear_obj_have_regularization:
             return self.reconstruction
 
diff --git a/autoarray/inversion/inversion/imaging/inversion_imaging_numba_util.py b/autoarray/inversion/inversion/imaging/inversion_imaging_numba_util.py
index 0eb95d26d..963e98e8e 100644
--- a/autoarray/inversion/inversion/imaging/inversion_imaging_numba_util.py
+++ b/autoarray/inversion/inversion/imaging/inversion_imaging_numba_util.py
@@ -857,3 +857,71 @@ def mapped_reconstructed_data_via_image_to_pix_unique_from(
             )
 
     return mapped_reconstructed_data
+
+
+@numba_util.jit()
+def relocated_grid_via_jit_from(grid, border_grid):
+    """
+    Relocate the coordinates of a grid to its border if they are outside the border, where the border is
+    defined as all pixels at the edge of the grid's mask (see *mask._border_1d_indexes*).
+
+    This is performed as follows:
+
+    1: Use the mean value of the grid's y and x coordinates to determine the origin of the grid.
+    2: Compute the radial distance of every grid coordinate from the origin.
+    3: For every coordinate, find its nearest pixel in the border.
+    4: Determine if it is outside the border, by comparing its radial distance from the origin to its paired
+    border pixel's radial distance.
+    5: If its radial distance is larger, use the ratio of radial distances to move the coordinate to the
+    border (if its inside the border, do nothing).
+
+    The method can be used on uniform or irregular grids, however for irregular grids the border of the
+    'image-plane' mask is used to define border pixels.
+
+    Parameters
+    ----------
+    grid
+        The grid (uniform or irregular) whose pixels are to be relocated to the border edge if outside it.
+    border_grid : Grid2D
+        The grid of border (y,x) coordinates.
+    """
+
+    grid_relocated = np.zeros(grid.shape)
+    grid_relocated[:, :] = grid[:, :]
+
+    border_origin = np.zeros(2)
+    border_origin[0] = np.mean(border_grid[:, 0])
+    border_origin[1] = np.mean(border_grid[:, 1])
+    border_grid_radii = np.sqrt(
+        np.add(
+            np.square(np.subtract(border_grid[:, 0], border_origin[0])),
+            np.square(np.subtract(border_grid[:, 1], border_origin[1])),
+        )
+    )
+    border_min_radii = np.min(border_grid_radii)
+
+    grid_radii = np.sqrt(
+        np.add(
+            np.square(np.subtract(grid[:, 0], border_origin[0])),
+            np.square(np.subtract(grid[:, 1], border_origin[1])),
+        )
+    )
+
+    for pixel_index in range(grid.shape[0]):
+        if grid_radii[pixel_index] > border_min_radii:
+            closest_pixel_index = np.argmin(
+                np.square(grid[pixel_index, 0] - border_grid[:, 0])
+                + np.square(grid[pixel_index, 1] - border_grid[:, 1])
+            )
+
+            move_factor = (
+                border_grid_radii[closest_pixel_index] / grid_radii[pixel_index]
+            )
+
+            if move_factor < 1.0:
+                grid_relocated[pixel_index, :] = (
+                    move_factor * (grid[pixel_index, :] - border_origin[:])
+                    + border_origin[:]
+                )
+
+    return grid_relocated
diff --git a/autoarray/inversion/inversion/imaging/mapping.py b/autoarray/inversion/inversion/imaging/mapping.py
index f72e27b04..7606e961e 100644
--- a/autoarray/inversion/inversion/imaging/mapping.py
+++ b/autoarray/inversion/inversion/imaging/mapping.py
@@ -1,6 +1,8 @@
 import numpy as np
 from typing import Dict, List, Optional, Union
 
+from autoconf import cached_property
+
 from autoarray.dataset.imaging.dataset import Imaging
 from autoarray.inversion.inversion.dataset_interface import DatasetInterface
 from autoarray.inversion.inversion.imaging.abstract import AbstractInversionImaging
@@ -89,7 +91,7 @@ def _data_vector_mapper(self) -> np.ndarray:
 
         return data_vector
 
-    @property
+    @cached_property
     def data_vector(self) -> np.ndarray:
         """
         The `data_vector` is a 1D vector whose values are solved for by the simultaneous linear equations constructed
@@ -156,7 +158,7 @@ def _curvature_matrix_mapper_diag(self) -> Optional[np.ndarray]:
 
         return curvature_matrix
 
-    @property
+    @cached_property
     def curvature_matrix(self):
         """
         The `curvature_matrix` is a 2D matrix which uses the mappings between the data and the linear objects to
diff --git a/autoarray/inversion/inversion/imaging/w_tilde.py b/autoarray/inversion/inversion/imaging/w_tilde.py
index 67d42544c..bc2fa8bf0 100644
--- a/autoarray/inversion/inversion/imaging/w_tilde.py
+++ b/autoarray/inversion/inversion/imaging/w_tilde.py
@@ -1,6 +1,8 @@
 import numpy as np
 from typing import Dict, List, Optional, Union
 
+from autoconf import cached_property
+
 from autoarray.dataset.imaging.dataset import Imaging
 from autoarray.dataset.imaging.w_tilde import WTildeImaging
 from autoarray.inversion.inversion.dataset_interface import DatasetInterface
@@ -64,9 +66,8 @@ def __init__(
 
         self.w_tilde = dataset.w_tilde
 
-    @property
+    @cached_property
     def w_tilde_data(self):
-
         return inversion_imaging_numba_util.w_tilde_data_imaging_from(
             image_native=np.array(self.data.native.array),
             noise_map_native=self.noise_map.native.array,
@@ -110,7 +111,7 @@ def _data_vector_mapper(self) -> np.ndarray:
 
         return data_vector
 
-    @property
+    @cached_property
     def data_vector(self) -> np.ndarray:
         """
         Returns the `data_vector`, a 1D vector whose values are solved for by the simultaneous linear equations
@@ -210,7 +211,7 @@ def _data_vector_func_list_and_mapper(self) -> np.ndarray:
 
         return data_vector
 
-    @property
+    @cached_property
     def curvature_matrix(self) -> np.ndarray:
         """
         Returns the `curvature_matrix`, a 2D matrix which uses the mappings between the data and the linear objects to
@@ -231,7 +232,6 @@ def curvature_matrix(self) -> np.ndarray:
         to ensure if we access it after computing the `curvature_reg_matrix` it is correctly recalculated in a new
         array of memory.
         """
-
         if self.has(cls=AbstractLinearObjFuncList):
             curvature_matrix = self._curvature_matrix_func_list_and_mapper
         elif self.total(cls=AbstractMapper) == 1:
diff --git a/autoarray/inversion/linear_obj/func_list.py b/autoarray/inversion/linear_obj/func_list.py
index ea7a583e4..b329a0c2f 100644
--- a/autoarray/inversion/linear_obj/func_list.py
+++ b/autoarray/inversion/linear_obj/func_list.py
@@ -1,6 +1,8 @@
 import numpy as np
 from typing import Optional
 
+from autoconf import cached_property
+
 from autoarray.inversion.linear_obj.linear_obj import LinearObj
 from autoarray.inversion.linear_obj.neighbors import Neighbors
 from autoarray.inversion.linear_obj.unique_mappings import UniqueMappings
@@ -44,7 +46,7 @@ def __init__(
 
         self.grid = grid
 
-    @property
+    @cached_property
     def neighbors(self) -> Neighbors:
         """
         An object describing how the different parameters in the linear object neighbor one another, which is used
@@ -76,7 +78,7 @@ def neighbors(self) -> Neighbors:
             arr=neighbors.astype("int"), sizes=neighbors_sizes.astype("int")
         )
 
-    @property
+    @cached_property
     def unique_mappings(self) -> UniqueMappings:
         """
         Returns the unique mappings of every unmasked data pixel's (e.g. `grid_slim`) sub-pixels (e.g. `grid_sub_slim`)
@@ -88,7 +90,6 @@ def unique_mappings(self) -> UniqueMappings:
 
         For a `LinearObjFuncList` every data pixel's group of sub-pixels maps directly to the linear function.
         """
-
         sub_size = np.max(self.grid.over_sample_size)
 
         # TODO : This shape slim is prob unreliable and needs to be divided by sub_size**2
diff --git a/autoarray/inversion/pixelization/border_relocator.py b/autoarray/inversion/pixelization/border_relocator.py
index 4856bfd10..d6dfe7545 100644
--- a/autoarray/inversion/pixelization/border_relocator.py
+++ b/autoarray/inversion/pixelization/border_relocator.py
@@ -2,6 +2,8 @@
 import numpy as np
 from typing import Tuple, Union
 
+from autoconf import cached_property
+
 from autoarray.mask.mask_2d import Mask2D
 from autoarray.structures.arrays.uniform_2d import Array2D
 from autoarray.structures.grids.uniform_2d import Grid2D
@@ -266,7 +268,9 @@ def relocated_grid_from(grid, border_grid, xp=np):
 
 
 class BorderRelocator:
-    def __init__(self, mask: Mask2D, sub_size: Union[int, Array2D]):
+    def __init__(
+        self, mask: Mask2D, sub_size: Union[int, Array2D], use_w_tilde: bool = False
+    ):
         """
         Relocates source plane coordinates that trace outside the mask’s border in the source-plane back onto the
         border.
@@ -323,6 +327,8 @@ def __init__(self, mask: Mask2D, sub_size: Union[int, Array2D]):
 
         self.sub_border_grid = sub_grid[self.sub_border_slim]
 
+        self.use_w_tilde = use_w_tilde
+
     def relocated_grid_from(self, grid: Grid2D, xp=np) -> Grid2D:
         """
         Relocate the coordinates of a grid to the border of this grid if they are outside the border, where the
@@ -350,21 +356,39 @@ def relocated_grid_from(self, grid: Grid2D, xp=np) -> Grid2D:
         if len(self.sub_border_grid) == 0:
             return grid
 
-        values = relocated_grid_from(
-            grid=grid.array, border_grid=grid.array[self.border_slim], xp=xp
-        )
+        if not self.use_w_tilde:
 
-        over_sampled = relocated_grid_from(
-            grid=grid.over_sampled.array,
-            border_grid=grid.over_sampled.array[self.sub_border_slim],
-            xp=xp,
-        )
+            values = relocated_grid_from(
+                grid=grid.array, border_grid=grid.array[self.border_slim], xp=xp
+            )
+
+            over_sampled = relocated_grid_from(
+                grid=grid.over_sampled.array,
+                border_grid=grid.over_sampled.array[self.sub_border_slim],
+                xp=xp,
+            )
+
+        else:
+
+            from autoarray.inversion.inversion.imaging import (
+                inversion_imaging_numba_util,
+            )
+
+            values = inversion_imaging_numba_util.relocated_grid_via_jit_from(
+                grid=grid.array, border_grid=grid[self.border_slim]
+            )
+
+            over_sampled = inversion_imaging_numba_util.relocated_grid_via_jit_from(
+                grid=grid.over_sampled.array,
+                border_grid=grid.over_sampled[self.sub_border_slim],
+            )
 
         return Grid2D(
             values=values,
             mask=grid.mask,
             over_sample_size=self.sub_size,
             over_sampled=over_sampled,
+            over_sampler=grid.over_sampler,
             xp=xp,
         )
 
@@ -384,9 +408,26 @@ def relocated_mesh_grid_from(
         if len(self.sub_border_grid) == 0:
             return mesh_grid
 
+        if not self.use_w_tilde:
+
+            relocated_grid = relocated_grid_from(
+                grid=mesh_grid.array,
+                border_grid=grid[self.border_slim],
+                xp=xp,
+            )
+
+        else:
+
+            from autoarray.inversion.inversion.imaging import (
+                inversion_imaging_numba_util,
+            )
+
+            relocated_grid = inversion_imaging_numba_util.relocated_grid_via_jit_from(
+                grid=mesh_grid.array,
+                border_grid=grid[self.sub_border_slim],
+            )
+
         return Grid2DIrregular(
-            values=relocated_grid_from(
-                grid=mesh_grid.array, border_grid=grid[self.sub_border_slim], xp=xp
-            ),
+            values=relocated_grid,
             xp=xp,
         )
diff --git a/autoarray/inversion/pixelization/mappers/abstract.py b/autoarray/inversion/pixelization/mappers/abstract.py
index 38800f4df..61dce22de 100644
--- a/autoarray/inversion/pixelization/mappers/abstract.py
+++ b/autoarray/inversion/pixelization/mappers/abstract.py
@@ -3,6 +3,7 @@
 from typing import List, Optional, Tuple
 
 from autoconf import conf
+from autoconf import cached_property
 
 from autoarray.inversion.linear_obj.linear_obj import LinearObj
 from autoarray.inversion.linear_obj.func_list import UniqueMappings
@@ -208,7 +209,7 @@ def sub_slim_indexes_for_pix_index(self) -> List[List]:
 
         return sub_slim_indexes_for_pix_index
 
-    @property
+    @cached_property
     def unique_mappings(self) -> UniqueMappings:
         """
         Returns the unique mappings of every unmasked data pixel's (e.g. `grid_slim`) sub-pixels (e.g. `grid_sub_slim`)
@@ -221,7 +222,6 @@ def unique_mappings(self) -> UniqueMappings:
         A full description of these mappings is given in the
         function `mapper_util.data_slim_to_pixelization_unique_from()`.
         """
-
         (
             data_to_pix_unique,
             data_weights,
@@ -245,7 +245,7 @@ def unique_mappings(self) -> UniqueMappings:
             pix_lengths=pix_lengths,
         )
 
-    @property
+    @cached_property
     def mapping_matrix(self) -> np.ndarray:
         """
         The `mapping_matrix` of a linear object describes the mappings between the observed data's data-points / pixels
diff --git a/autoarray/inversion/pixelization/mappers/delaunay.py b/autoarray/inversion/pixelization/mappers/delaunay.py
index 3be7dd97b..1fe990bee 100644
--- a/autoarray/inversion/pixelization/mappers/delaunay.py
+++ b/autoarray/inversion/pixelization/mappers/delaunay.py
@@ -1,5 +1,7 @@
 import numpy as np
 
+from autoconf import cached_property
+
 from autoarray.inversion.pixelization.mappers.abstract import AbstractMapper
 from autoarray.inversion.pixelization.mappers.abstract import PixSubWeights
 
@@ -59,7 +61,7 @@ class MapperDelaunay(AbstractMapper):
     def delaunay(self):
         return self.source_plane_mesh_grid.delaunay
 
-    @property
+    @cached_property
     def pix_sub_weights(self) -> PixSubWeights:
         """
         Computes the following three quantities describing the mappings between of every sub-pixel in the masked data
diff --git a/autoarray/inversion/pixelization/mappers/rectangular.py b/autoarray/inversion/pixelization/mappers/rectangular.py
index f8c227c4d..8662ea38e 100644
--- a/autoarray/inversion/pixelization/mappers/rectangular.py
+++ b/autoarray/inversion/pixelization/mappers/rectangular.py
@@ -1,5 +1,7 @@
 from typing import Tuple
 
+from autoconf import cached_property
+
 from autoarray.inversion.pixelization.mappers.abstract import AbstractMapper
 from autoarray.inversion.pixelization.mappers.abstract import PixSubWeights
 
@@ -58,7 +60,7 @@ class MapperRectangular(AbstractMapper):
     def shape_native(self) -> Tuple[int, ...]:
         return self.source_plane_mesh_grid.shape_native
 
-    @property
+    @cached_property
     def pix_sub_weights(self) -> PixSubWeights:
         """
         Computes the following three quantities describing the mappings between of every sub-pixel in the masked data
diff --git a/autoarray/inversion/pixelization/mappers/voronoi.py b/autoarray/inversion/pixelization/mappers/voronoi.py
index cef9bf36c..62b03afc2 100644
--- a/autoarray/inversion/pixelization/mappers/voronoi.py
+++ b/autoarray/inversion/pixelization/mappers/voronoi.py
@@ -1,6 +1,9 @@
 import numpy as np
 from typing import Optional, Tuple
 
+
+from autoconf import cached_property
+
 from autoarray.inversion.pixelization.mappers.abstract import AbstractMapper
 from autoarray.inversion.pixelization.mappers.abstract import PixSubWeights
 from autoarray.structures.arrays.uniform_2d import Array2D
@@ -80,7 +83,7 @@ def pix_sub_weights_split_cross(self) -> PixSubWeights:
 
         return PixSubWeights(mappings=mappings, sizes=sizes, weights=weights)
 
-    @property
+    @cached_property
     def pix_sub_weights(self) -> PixSubWeights:
         """
         Computes the following three quantities describing the mappings between of every sub-pixel in the masked data
diff --git a/autoarray/inversion/pixelization/mesh/abstract.py b/autoarray/inversion/pixelization/mesh/abstract.py
index 2dceeb590..f277e49a1 100644
--- a/autoarray/inversion/pixelization/mesh/abstract.py
+++ b/autoarray/inversion/pixelization/mesh/abstract.py
@@ -47,6 +47,7 @@ def relocated_grid_from(
             mask=source_plane_data_grid.mask,
             over_sample_size=source_plane_data_grid.over_sampler.sub_size,
             over_sampled=source_plane_data_grid.over_sampled.array,
+            over_sampler=source_plane_data_grid.over_sampler,
             xp=xp,
         )
 
diff --git a/autoarray/operators/over_sampling/over_sampler.py b/autoarray/operators/over_sampling/over_sampler.py
index e9942547e..a699eae60 100644
--- a/autoarray/operators/over_sampling/over_sampler.py
+++ b/autoarray/operators/over_sampling/over_sampler.py
@@ -3,6 +3,7 @@
 from typing import Union
 
 from autoconf import conf
+from autoconf import cached_property
 
 from autoarray.mask.mask_2d import Mask2D
 from autoarray.structures.arrays.uniform_2d import Array2D
@@ -330,7 +331,7 @@ def slim_for_sub_slim(self) -> np.ndarray:
             mask_2d=np.array(self.mask), sub_size=self.sub_size.array
         ).astype("int")
 
-    @property
+    @cached_property
     def uniform_over_sampled(self):
         """
         For a sub-grid, every unmasked pixel of its 2D mask with shape (total_y_pixels, total_x_pixels) is divided into
@@ -345,6 +346,7 @@ def uniform_over_sampled(self):
         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.
         """
+
         from autoarray.structures.grids.irregular_2d import Grid2DIrregular
 
         grid = over_sample_util.grid_2d_slim_over_sampled_via_mask_from(
diff --git a/autoarray/structures/decorators/to_grid.py b/autoarray/structures/decorators/to_grid.py
index e2d36fda1..acad17bbe 100644
--- a/autoarray/structures/decorators/to_grid.py
+++ b/autoarray/structures/decorators/to_grid.py
@@ -26,17 +26,25 @@ def via_grid_2d(self, result) -> Union[Grid2D, List[Grid2D]]:
             except AttributeError:
                 over_sampled = None
 
+            try:
+                over_sampler = result.over_sampler
+            except AttributeError:
+                over_sampler = None
+
             return Grid2D(
                 values=result,
                 mask=self.mask,
                 over_sample_size=self.over_sample_size,
                 over_sampled=over_sampled,
+                over_sampler=over_sampler,
             )
 
         try:
             grid_over_sampled_list = [res.over_sampled for res in result]
+            grid_over_sampler_list = [res.over_sampler for res in result]
         except AttributeError:
             grid_over_sampled_list = [None] * len(result)
+            grid_over_sampler_list = [None] * len(result)
 
         return [
             Grid2D(
@@ -44,8 +52,11 @@ def via_grid_2d(self, result) -> Union[Grid2D, List[Grid2D]]:
                 mask=self.mask,
                 over_sample_size=self.over_sample_size,
                 over_sampled=over_sampled,
+                over_sampler=over_sampler,
+            )
+            for res, over_sampled, over_sampler in zip(
+                result, grid_over_sampled_list, grid_over_sampler_list
             )
-            for res, over_sampled in zip(result, grid_over_sampled_list)
         ]
 
     def via_grid_2d_irr(self, result) -> Union[Grid2DIrregular, List[Grid2DIrregular]]:
diff --git a/autoarray/structures/grids/uniform_2d.py b/autoarray/structures/grids/uniform_2d.py
index 6292549b3..4dc73b498 100644
--- a/autoarray/structures/grids/uniform_2d.py
+++ b/autoarray/structures/grids/uniform_2d.py
@@ -27,6 +27,7 @@ def __init__(
         store_native: bool = False,
         over_sample_size: Union[int, Array2D] = 4,
         over_sampled: Optional[Grid2D] = None,
+        over_sampler=None,
         xp=np,
         *args,
         **kwargs,
@@ -169,15 +170,36 @@ def __init__(
 
         grid_2d_util.check_grid_2d(grid_2d=values)
 
-        over_sample_size = over_sample_util.over_sample_size_convert_to_array_2d_from(
-            over_sample_size=over_sample_size, mask=mask
+        self._over_sample_size = over_sample_size
+
+        self._over_sampler = over_sampler
+        self._over_sampled = over_sampled
+
+    @property
+    def over_sample_size(self):
+
+        if isinstance(self._over_sample_size, Array2D):
+            return self._over_sample_size
+
+        self._over_sample_size = (
+            over_sample_util.over_sample_size_convert_to_array_2d_from(
+                over_sample_size=self._over_sample_size, mask=self.mask
+            )
         )
 
+        return self._over_sample_size
+
+    @property
+    def over_sampler(self):
+
+        if self._over_sampler is not None:
+            return self._over_sampler
+
         from autoarray.operators.over_sampling.over_sampler import OverSampler
 
-        self.over_sampler = OverSampler(sub_size=over_sample_size, mask=mask)
+        self._over_sampler = OverSampler(sub_size=self.over_sample_size, mask=self.mask)
 
-        self._over_sampled = over_sampled
+        return self._over_sampler
 
     @property
     def over_sampled(self):
@@ -188,7 +210,7 @@ def over_sampled(self):
         over_sampled = over_sample_util.grid_2d_slim_over_sampled_via_mask_from(
             mask_2d=np.array(self.mask),
             pixel_scales=self.mask.pixel_scales,
-            sub_size=self.over_sampler.sub_size.array.astype("int"),
+            sub_size=self.over_sample_size.array.astype("int"),
             origin=self.mask.origin,
         )
 
@@ -700,12 +722,9 @@ def subtracted_from(
             mask=mask,
             over_sample_size=self.over_sample_size,
             over_sampled=self.over_sampled - xp.array(offset),
+            over_sampler=self.over_sampler,
         )
 
-    @property
-    def over_sample_size(self):
-        return self.over_sampler.sub_size
-
     @property
     def slim(self) -> "Grid2D":
         """