# Using the GridR's core grid resampling API This guide shows how to use the core `array_grid_resampling` GridR's feature. It covers the following topics: - Performing basic grid-based resampling - Working with spatial windows - Controlling output grid resolution - Handling nodata values and generating output validity masks - Moving the input data origin - Using masks to control the validity of input grids and input data ## Setting things up First we have to import the `array_grid_resampling` method. We also import the well known mandrill raster image. ```python from gridr.core.grid.grid_resampling import array_grid_resampling from gridr.misc.mandrill import mandrill ``` Please note the mandrill image has 3 color's channels. Here we limit the image to the first channel and convert its type to 64 bits floating points. ```python image = mandrill[0] grid_dtype = np.float64 data_in_dtype = np.float64 data_out_dtype = np.float64 array_in = image.astype(data_in_dtype) ``` ## Identity grid transform We start by defining the idendity transform grid. That's not really impressive but we will go further to illustrate some variants (windowing, zooming, ...) Let's first create the idendity transformation grid. ```python # create identity grid if image.ndim == 2: x = np.arange(0, image.shape[0], dtype=grid_dtype) y = np.arange(0, image.shape[1], dtype=grid_dtype) xx, yy = np.meshgrid(x, y) ``` ### Simple idendity transform The `array_grid_resampling` method provides the `array_out` argument in order to use a preallocated C-contiguous array as output. It is quite convenient when implementing a tiling chain that uses that method. For now, we choose to keep it simple : by setting `array_out` to `None` we let the method allocate the output buffer for us. We have defined the rows and columns grids to be at full resolution, ie. there is one coordinate in the grids for each pixel in the input image. To stay at that resolution we set `grid_resolution` to 1 in both directions. The `array_grid_resampling` function returns a tuple of two NumPy arrays. For now, we will focus on the first one — it contains the allocated output buffer, which holds the result of the interpolation. For the interpolation method, we will employ the `cubic` identifier, which designates the Optimized Bicubic interpolation algorithm. ```python # Lets call the grid resampling array_out_id, _ = array_grid_resampling( interp="cubic", array_in=image.astype(data_in_dtype), grid_row=yy, grid_col=xx, grid_resolution=(1,1), array_out=None, win=None, array_in_mask=None, grid_mask=None, array_out_mask=None, nodata_out=None, ) ``` ![001.png](grid_resampling_001_work_files/001.png) In the following test we compare input with output. It must succeed ! ```python np.testing.assert_almost_equal(array_in, array_out_id, decimal=10) ``` ### Apply a window The `array_grid_resampling` method provides the `win` parameter in order to set a window to limit the region to compute. The indices used in the `win` definition refers to indices considering the full resolution output. Below, we first define a window centered on the mandrill's left eye. Since we are considering an identity grid, the grid coordinates align with the input array's image coordinates ```python win_center = np.array((58, 175)) # left eye win_shape = np.array((100, 100)) win = np.array((win_center - win_shape//2, win_center - win_shape//2 + win_shape - 1)).T ``` ![002.png](grid_resampling_001_work_files/002.png) ```python array_out_id_win_1_1, _ = array_grid_resampling( interp="cubic", array_in=array_in, grid_row=yy, grid_col=xx, grid_resolution=(1,1), array_out=None, win=win, array_in_mask=None, grid_mask=None, array_out_mask=None, nodata_out=None, ) ``` ![003.png](grid_resampling_001_work_files/003.png) ### Playing with resolution to apply a zoom The `array_grid_resampling` method uses the `grid_resolution` parameter to insert additional nodes between existing grid nodes. This insertion is achieved through bilinear interpolation, with the oversampling factors for both rows and columns specified in the grid_resolution parameter. While the primary purpose of the `grid_resolution` is to control computation time and grid storage, in this case, we use it to perform a zoom operation, increasing the row resolution by a factor of 2 and the column resolution by a factor of 5. At the same time, we want to maintain the same input region as before. To achieve this, the window definition must be adjusted to accommodate the new full-resolution output. ```python # Set resolution to the target zoom factors resolution = (2, 3) # The full resolution grid has (Nrow-1)*resolution_row + 1 rows (apply same logic on column) full_output_shape = (np.asarray(array_in.shape) - 1)*resolution + 1 # setting output window to only consider a region of interest # here we construct a window centered on the left eye in the full resolution output grid win_center = np.array((58, 175)) * np.array(resolution) win_shape = np.array((100, 100)) * np.array(resolution) win = np.array((win_center - win_shape//2, win_center - win_shape//2 + win_shape - 1)).T # Lets call the grid resampling array_out_id_zoom_2_5, _ = array_grid_resampling( interp="cubic", array_in=array_in, grid_row=yy, grid_col=xx, grid_resolution=resolution, array_out=None, win=win, array_in_mask=None, grid_mask=None, array_out_mask=None, nodata_out=None, ) ``` ![004.png](grid_resampling_001_work_files/004.png) In order to illustrate the change of interpolation method, we will use here the `nearest` identifier in order to adress a nearest neighbor interpolation. ```python # Lets call the grid resampling array_out_id_zoom_2_5_nearest, _ = array_grid_resampling( interp="nearest", array_in=array_in, grid_row=yy, grid_col=xx, grid_resolution=resolution, array_out=None, win=win, array_in_mask=None, grid_mask=None, array_out_mask=None, nodata_out=None, ) # Let's show it plot_im(array_out_id_zoom_2_5_nearest, None, prefix='004_nearest') ``` ![004_nearest.png](grid_resampling_001_work_files/004_nearest.png) ## Geometry transformation examples In that part we provide some examples of simple non identity transformation. ### Translation transform Let's begin with a simple translation on the column coordinates by shifting the image to the left. ```python # First create identity grid if image.ndim == 2: x = np.arange(0, image.shape[0], dtype=grid_dtype) y = np.arange(0, image.shape[1], dtype=grid_dtype) xx, yy = np.meshgrid(x, y) xx += 50.5 # Lets call the grid resampling array_out_translation_left, _ = array_grid_resampling( interp="cubic", array_in=array_in, grid_row=yy, grid_col=xx, grid_resolution=(1,1), array_out=None, win=None, array_in_mask=None, grid_mask=None, array_out_mask=None, nodata_out=None, ) ``` ![005.png](grid_resampling_001_work_files/005.png) ```python array_out_translation_left[:, -53:] ``` array([[124.875 , nan, nan, ..., nan, nan, nan], [130. , nan, nan, ..., nan, nan, nan], [102.8125, nan, nan, ..., nan, nan, nan], ..., [102.5625, nan, nan, ..., nan, nan, nan], [ 86.5625, nan, nan, ..., nan, nan, nan], [ 3.5625, nan, nan, ..., nan, nan, nan]], shape=(512, 53)) As you can see, the right strip border is filled with `NaN` values. It is due to the fact that the grid target coordinates went outside of the image's definition domain. Please note that the filling with `NaN` is not the default behaviour and it occured because we set the `nodata_out` parameter to `None`. Let's replay it by setting `nodata_out` to its default value, ie 0. As you will see, the right strip border is now filled with zeros. We also take this opportunity to introduce the use of the `array_mask_out` parameter and examine the second array returned by the function. We will verify that the mask correctly matches the `nodata` region of the resampled image. Here we set the `array_out_mask` parameter to `True` ; that's a convenient way to let the `array_grid_resampling` method allocate the right sized `unsigned byte` (ie. uint8) array and return it. Please note that we can also directly pass a preallocated array similarly to the `array_out` parameter. ```python array_out_translation_left, array_out_translation_left_mask = array_grid_resampling( interp="cubic", array_in=array_in, grid_row=yy, grid_col=xx, grid_resolution=(1,1), array_out=None, win=None, array_in_mask=None, grid_mask=None, array_out_mask=True, nodata_out=0., ) array_out_translation_left[:, -53:] ``` array([[124.875 , 0. , 0. , ..., 0. , 0. , 0. ], [130. , 0. , 0. , ..., 0. , 0. , 0. ], [102.8125, 0. , 0. , ..., 0. , 0. , 0. ], ..., [102.5625, 0. , 0. , ..., 0. , 0. , 0. ], [ 86.5625, 0. , 0. , ..., 0. , 0. , 0. ], [ 3.5625, 0. , 0. , ..., 0. , 0. , 0. ]], shape=(512, 53)) Let's display the output validity mask ```python array_out_translation_left_mask[:, -53:] ``` array([[1, 0, 0, ..., 0, 0, 0], [1, 0, 0, ..., 0, 0, 0], [1, 0, 0, ..., 0, 0, 0], ..., [1, 0, 0, ..., 0, 0, 0], [1, 0, 0, ..., 0, 0, 0], [1, 0, 0, ..., 0, 0, 0]], shape=(512, 53), dtype=uint8) ![005-mask-out.png](grid_resampling_001_work_files/005-mask-out.png) ### Rotation transform Here we define a grid in order to apply a 45 degree rotation. In order to have the full image rotated please note that we have to extend the grid size. ```python center = np.asarray(array_in.shape) // 2 theta = np.pi/4. # we know the image is square - here we extend the output grid in order to cover # the full image ext = image.shape[0] * (2**0.5 - 1) // 2 # First create identity grid if image.ndim == 2: x = np.arange(0, image.shape[0] + 2 * ext, dtype=grid_dtype) - ext y = np.arange(0, image.shape[1] + 2 * ext, dtype=grid_dtype) - ext xx, yy = np.meshgrid(x , y) # Apply the rotation xrot = (xx - center[1]) * np.cos(theta) - (yy - center[0]) * np.sin(theta) + center[1] yrot = (xx - center[1]) * np.sin(theta) + (yy - center[0]) * np.cos(theta) + center[0] ``` As a demonstration purpose we here preallocate the array_out buffer. ```python # We have to init the out buffer with the right shape array_out = np.zeros(xx.shape, dtype=data_out_dtype) # Lets call the grid resampling _, mask_out = array_grid_resampling( interp="cubic", array_in=array_in, grid_row=yrot, grid_col=xrot, grid_resolution=(1,1), array_out=array_out, win=None, array_in_mask=None, grid_mask=None, array_out_mask=True, nodata_out=None, ) ``` ![006.png](grid_resampling_001_work_files/006.png) ### Specifying an output target window In the previous step we opt to pre-allocate the output buffer for the image. What if the requested production window is smaller that the `array_out` buffer ? ```python # Set a window to focus on the left eye (we get the coordinates from the previous resampled image) win = np.asarray([[250, 320], [120, 220]]) # We have to init the out buffer with an sufficient shape array_out = np.zeros((200, 200), dtype=data_out_dtype) # Lets call the grid resampling _, mask_out = array_grid_resampling( interp="cubic", array_in=array_in, grid_row=yrot, grid_col=xrot, grid_resolution=(1,1), array_out=array_out, win=win, array_in_mask=None, grid_mask=None, array_out_mask=True, nodata_out=None, ) ``` ![007.png](grid_resampling_001_work_files/007.png) It’s worth noting that even without passing a preallocated output validity mask, the internally generated one will still match the shape of the `array_out` buffer. As shown above, the resampled window is written starting at the origin of `array_out` (i.e., at [0, 0]). In some cases - for example, when working with tiles that share a common output buffer - it may be necessary to write the result to a specific region within `array_out`. This is precisely the purpose of the `array_out_win` parameter. ```python # Set a window to focus on the left eye (we get the coordinates from the previous resampled image) win = np.asarray([[250, 320], [120, 220]]) # We have to init the out buffer with an sufficient shape # Here we will also demonstrate that the resampling only modify the target window # For that purpose we do not initialize with zeros but with arange array_out = np.arange(200*200, dtype=data_out_dtype).reshape((200, 200)) % 255 # We also preallocate the validity mask for the same purpose mask_out = np.arange(array_out.size, dtype=np.uint8).reshape(array_out.shape) % 10 # Let's define where we want to place the resampled window array_out_win_origin = np.asarray([25, 40]) array_out_win = np.asarray([ [array_out_win_origin[0], array_out_win_origin[0] + win[0][1] - win[0][0]], [array_out_win_origin[1], array_out_win_origin[1] + win[1][1] - win[1][0]]]) # Lets call the grid resampling _, _ = array_grid_resampling( interp="cubic", array_in=array_in, grid_row=yrot, grid_col=xrot, grid_resolution=(1,1), array_out=array_out, array_out_win=array_out_win, # <= new parameter win=win, array_in_mask=None, grid_mask=None, array_out_mask=mask_out, nodata_out=None, ) ``` ![008.png](grid_resampling_001_work_files/008.png) Note that this example also demonstrates that operating on a specific window only affects the corresponding region in the output array, leaving all other areas unchanged. ## Adjusting the Origin of Input Coordinates The `array_grid_resampling` function includes the parameter `array_in_origin` (default: (0., 0.)) that allows shifting the origin of the input coordinate system. From an implementation standpoint, this shift is applied directly to the computed target grid coordinates. This feature is useful in several scenarios: - Handling different grid conventions without modifying the grid itself. As a reminder, GridR assumes that integer coordinates refer to pixel centers. In other systems, a pixel may be treated as an area, with its center located relatively at (0.5, 0.5). - Working with subregions of a larger input array. In this case, the shift can represent the position of the subregion within the full array’s coordinate system. As an example we will take a previous example using an idendity transformation. ```python # create identity grid if image.ndim == 2: x = np.arange(0, image.shape[0], dtype=grid_dtype) y = np.arange(0, image.shape[1], dtype=grid_dtype) xx, yy = np.meshgrid(x, y) win_center = np.array((58, 175)) # left eye win_shape = np.array((100, 100)) win = np.array((win_center - win_shape//2, win_center - win_shape//2 + win_shape - 1)).T # Call it without specifying array_in_origin array_out_origin_0_0, _ = array_grid_resampling( interp="cubic", array_in=array_in, grid_row=yy, grid_col=xx, grid_resolution=(1,1), array_out=None, win=win, array_in_mask=None, grid_mask=None, array_out_mask=None, nodata_out=None, ) # Call it with specifying array_in_origin array_out_origin_30_m10, _ = array_grid_resampling( interp="cubic", array_in=array_in, array_in_origin=(30, -10), # <= shift +30 on rows and -10 on columns grid_row=yy, grid_col=xx, grid_resolution=(1,1), array_out=None, win=win, array_in_mask=None, grid_mask=None, array_out_mask=None, nodata_out=None, ) ``` ![009-array-origin.png](grid_resampling_001_work_files/009-array-origin.png) As you can observe, the coordinate shift is reflected in the bottom image. ## Working with masks In this part we will see how to use masking. The `array_grid_resampling` method provides several kind of masking for both the grid and the input array. ### General convention GridR resampling functionalities uses the following convention for mask : - the mask are considered unsigned 8 bits integers taking value in {0, 1}. - a value of 0 is considered as non valid - a value of 1 is considered as valid ### Grid masking #### Grid masking with a raster mask Here we use an 8 bit raster having the same shape of the grids in order to devalidate points in the grid. To keep it simple we will demonstrate using an identity transform. ```python # First create identity grid if image.ndim == 2: x = np.arange(0, image.shape[0], dtype=grid_dtype) y = np.arange(0, image.shape[1], dtype=grid_dtype) xx, yy = np.meshgrid(x, y) ``` We define the mask as following : * the first 10 lines and last 20 lines are masked * the first 30 columns and last 40 columns are masked * a 80x80 window positionned at the grid center are masked Please note the mask applies to the grid and not to the input image. By convention, the mask is considered binary. A value of 0 corresponds to a masked cell, a value of 1 corresponds to a valid cell. ```python # setting a mask to devalidate margins and center grid_mask_valid_value = 1 grid_mask_invalid_value = 0 # the final multiplication here is not necessary but we keep it to show the mechanic grid_mask = np.ones(xx.shape, dtype=np.uint8) * grid_mask_valid_value # masked pixels are positiv ones grid_mask[0:10] = grid_mask_invalid_value grid_mask[-20:] = grid_mask_invalid_value grid_mask[:,0:30] = grid_mask_invalid_value grid_mask[:,-40:] = grid_mask_invalid_value # devalidate center grid_mask[image.shape[0]//2-40:image.shape[0]//2+40, image.shape[1]//2-40:image.shape[1]//2+40] = grid_mask_invalid_value ``` The `nodata_out` is used to set the output value affected by the masking operation. ```python # Lets call the grid resampling array_out, mask_out = array_grid_resampling( interp="cubic", array_in=array_in, grid_row=yy, grid_col=xx, grid_resolution=(1,1), array_out=None, win=None, array_in_mask=None, grid_mask=grid_mask, grid_mask_valid_value=grid_mask_valid_value, array_out_mask=True, nodata_out=0, ) ``` ![010.png](grid_resampling_001_work_files/010.png) As in the previous example, we pass True to the function to obtain the validity mask output. We then check that it aligns correctly with the grid mask. #### Masking with an inplaced grid nodata value In this case, we avoid providing an additional mask, primarily to limit memory usage. The `array_grid_resampling` method allows specifying a `grid_nodata` parameter, which is used to invalidate grid cells. Cells in the grid whose values match this nodata value (within a fixed absolute tolerance of 1e-5) are treated as masked. Providing a non-`None` value for `grid_nodata` activates the masking mechanism. Note that it is not necessary to assign the nodata value to both grid coordinate arrays (row and column); setting it in either one is sufficient for the corresponding cell to be considered masked. ```python # Set a nodata value grid_nodata = -999. # Apply it to a window of the grid yy[50:200,100:350] = grid_nodata ``` ```python # Lets call the grid resampling # Spoiler alert: this will fail ! try: array_out, _ = array_grid_resampling( interp="cubic", array_in=array_in, grid_row=yy, grid_col=xx, grid_resolution=(1,1), array_out=None, win=None, array_in_mask=None, grid_nodata=grid_nodata, grid_mask=grid_mask, grid_mask_valid_value=grid_mask_valid_value, array_out_mask=None, nodata_out=0, ) except ValueError as e: display(e) ``` ValueError('Only one of `grid_mask` or `grid_nodata` may be provided, not both.') The call failed, and the ValueError exception provides the reason: both masking modes, grid_mask and grid_nodata, cannot be used simultaneously. ```python # Lets remove the grid_mask parameter's definition array_out, _ = array_grid_resampling( interp="cubic", array_in=array_in, grid_row=yy, grid_col=xx, grid_resolution=(1,1), array_out=None, win=None, array_in_mask=None, grid_nodata=grid_nodata, grid_mask=None, grid_mask_valid_value=None, array_out_mask=None, nodata_out=0, ) ``` ![011.png](grid_resampling_001_work_files/011.png) #### Effect of Grid Resolution on Grid Masking Both the raster mask and the invalid_value grid mask remain consistent with the sampling of the grid. When the resolution differs from 1 in a given direction, interpolation must occur within each grid cell — that is, within the mesh defined by four neighboring grid nodes. In such cases, if at least one non-zero weighted node of the mesh is masked, the interpolated result is marked as invalid. This results in nodata values in the resampled image and corresponding invalid entries (value 0) in the output mask. Let's illustrate this behavior by zooming in on a small region. ```python # create identity grid on the eye if image.ndim == 2: x = np.arange(130, 220, dtype=grid_dtype) y = np.arange(30, 100, dtype=grid_dtype) xx, yy = np.meshgrid(x, y) resolution = np.asarray((2,3)) # lets define the mask using grid_nodata grid_nodata = -99999 yy[30, 35] = grid_nodata array_out, mask_out = array_grid_resampling( interp="cubic", array_in=array_in, grid_row=yy, grid_col=xx, grid_resolution=tuple(resolution), array_out=None, win=None, array_in_mask=None, grid_nodata=grid_nodata, array_out_mask=True, nodata_out=0, ) ``` ![012.png](grid_resampling_001_work_files/012.png) To invalidate a specific grid position, we've configured the input (30, 35). This grid coordinate directly corresponds to the (60, 105) output position in the output coordinate system. ```python print(mask_out[60-2:60+3,105-3:105+4]) print(f"Number of masked data : {np.where(mask_out==0)[0].shape}") ``` [[1 1 1 1 1 1 1] [1 0 0 0 0 0 1] [1 0 0 0 0 0 1] [1 0 0 0 0 0 1] [1 1 1 1 1 1 1]] Number of masked data : (15,) The masked area, as illustrated, exceeds a single pixel due to a resolution factor of 2 for rows and 3 for columns. The displayed array is centered on the position of the invalid grid node. Notably, this invalidation propagates to its immediate neighborhood, but exclusively for interpolated values where its weight is non-zero. For clearer comprehension, in the schema below, an `x` denotes an invalid grid node, an `o` represents other grid nodes, and both `-` and `*` signify oversampled positions. It's crucial to note that only the `*` positions utilize the `x` (invalid grid node) in their bilinear interpolation calculations. o - - o - - o - * * * * * - o * * x * * o - * * * * * - o - - o - - o ### Input Array masking #### Array masking with a raster mask Here we use an 8 bit raster having the same shape of the input array in order to devalidate points in the input array. Basically if input data is devalidate, interpolation that involves those data must be considered as not valid resulting in output valued to `nodata_out` and the optional `array_out_mask` set to 0 instead of 1. To keep it simple we will demonstrate using a simple zoom transformation transform. ```python array_in_mask = np.ones(array_in.shape, dtype=np.uint8) masked_pos = (60, 160) array_in_mask[*masked_pos] = 0 # In order to check the good interpolation we will also set the array_in to hold an aberrant value save_array_in = array_in[*masked_pos] # Create a grid centered on the masked position y = masked_pos[0] + np.linspace(-3, 3, 100) x = masked_pos[1] + np.linspace(-3, 3, 100) xx, yy = np.meshgrid(x, y) # First compute a array_out_womask_truth, _ = array_grid_resampling( interp="cubic", array_in=array_in, grid_row=yy, grid_col=xx, grid_resolution=(1,1), array_out=None, nodata_out=0, ) # Second compute without array masking to watch the effect of the non masked aberrant value array_in[*masked_pos] = 1e7 array_out_womask, _ = array_grid_resampling( interp="cubic", array_in=array_in, grid_row=yy, grid_col=xx, grid_resolution=(1,1), array_out=None, array_out_mask=True, nodata_out=0, ) aberrant_radiation = np.zeros(array_out_womask.shape) aberrant_radiation[np.where(array_out_womask - array_out_womask_truth != 0)] = 1 # Recompute with the mask array_out, mask_out = array_grid_resampling( interp="cubic", array_in=array_in, grid_row=yy, grid_col=xx, grid_resolution=(1,1), array_in_mask=array_in_mask, array_out=None, array_out_mask=True, nodata_out=150, # Here we set nodata_out to 150 in order to have some texture for the display ) # restore array_in array_in[*masked_pos] = save_array_in ``` ![013-array-masked-data.png](grid_resampling_001_work_files/013-array-masked-data.png) We have several plots to analyze here. The first plot reveals what happens when the anomalous value is not masked. While the 8-bit display dynamic range prevents us from fully appreciating the valid interpolation, this plot effectively showcases the influence of an anomalous value (exaggerated for illustrative purposes) on the interpolation process. The second plot, a binary representation, indicates where the anomalous value "spreads" or "irradiates." Interestingly, some unaffected lines can be seen traversing the "irradiated" zone. The third plot demonstrates the outcome when the anomalous value is successfully masked. We can see a clear correspondence between the no-data pixels and the "irradiated" area. Lastly, the final plot presents the output validity mask, confirming its consistency with our earlier findings. ```python # min and max indices of masked data upper_left = np.min(np.where(mask_out==0), axis=1) lower_right = np.max(np.where(mask_out==0), axis=1) print(f'upper left corner of masked data : {upper_left}') print(f'lower right corner of masked data : {lower_right}') # corresponding grid coordinates print(f'masked data upper left corner grid coordinates : {yy[*upper_left], xx[*upper_left]}') print(f'masked data lower right corner grid coordinates : {yy[*lower_right], xx[*lower_right]}') ``` upper left corner of masked data : [17 17] lower right corner of masked data : [82 82] masked data upper left corner grid coordinates : (np.float64(58.03030303030303), np.float64(158.03030303030303)) masked data lower right corner grid coordinates : (np.float64(61.96969696969697), np.float64(161.96969696969697)) The optimized bicubic interpolator uses a spatial support of ±2 around the target coordinates. This directly explains the masked area's boundaries: row indices [58,62] and column indices [158,162]. Now, let's discuss the unmasked lines that cross the invalid area. To do this, we will plot the grid coordinates, differentiating between the valid and invalid ones. ![013-mask-out-grid-scatter.png](grid_resampling_001_work_files/013-mask-out-grid-scatter.png) We've used a distinct style in the plot above to highlight coordinates that fall on exact integer values. This specific scenario, where we have whole integer coordinates within our grid, serves to prove that the `GridR` implementation aligns correctly with the "irradiated" area. It also shows that it avoids inadvertently masking pixels that should remain unmasked. What makes these integer coordinate positions unique? For these particular points, the subsequent interpolation process can utilize a directly available value in at least one direction, meaning no actual interpolation is required along that axis. In these instances, the neighboring points that would normally contribute to the interpolation—especially any anomalous ones—are effectively bypassed. Their associated coefficients are set to zero, which then yields a perfectly valid interpolated value.