Floating point error in feature detection breaks interpn

[travis-j-hahn]

It looks like there is a chance during the feature detection step for the identification of cells on boundary edges which result in hdim1 or hdim2 values slightly outside of the maximum boundary. I.e. the feature detection can create a cell with ix value of 218.00000003 when the maximum index is 218, and this throws an error in tobac/general.py:107 with coordinate_points = interpn(points, values, xi) raising ValueError: One of the requested xi is out of bounds in dimension 0

Full Error:

Traceback (most recent call last):

  Cell In[282], line 1
    wrf_features = wrf_tobac_feature_id(wrf_cube, 'IC', CONFIG)

  File ~/CoMET/CoMET/wrf_tobac.py:40 in wrf_tobac_feature_id
    wrf_radar_features = tobac.feature_detection.feature_detection_multithreshold(cube, dxy=dxy, **CONFIG['wrf']['tobac']['feature_id'])

  File ~/.conda/envs/CoMET/lib/python3.11/site-packages/tobac/utils/decorators.py:271 in wrapper
    output = func(*args, **kwargs)

  File ~/.conda/envs/CoMET/lib/python3.11/site-packages/tobac/feature_detection.py:1409 in feature_detection_multithreshold
    features = add_coordinates(features, field_in)

  File ~/.conda/envs/CoMET/lib/python3.11/site-packages/tobac/utils/general.py:109 in add_coordinates
    coordinate_points = interpn(points, values, xi)

  File ~/.conda/envs/CoMET/lib/python3.11/site-packages/scipy/interpolate/_rgi.py:740 in interpn
    raise ValueError("One of the requested xi is out of bounds "

ValueError: One of the requested xi is out of bounds in dimension 0

This can be fixed by using np.clip to ensure there is no overflow:

elif variable_cube.coord(coord).ndim == 2:
            if variable_cube.coord_dims(coord) == (hdim_1, hdim_2):
                points = (dimvec_1, dimvec_2)
                values = variable_cube.coord(coord).points
                temp_hdim1 = np.clip(t['hdim_1'].values, 0, points[0].max())
                temp_hdim2 = np.clip(t['hdim_2'].values, 0, points[1].max())
                xi = np.column_stack((temp_hdim1, temp_hdim2))
                coordinate_points = interpn(points, values, xi)