How to use the pymedphys._imports.numpy.max function in pymedphys

def create_dose_function(net_od, dose):
    net_od = np.array(net_od, copy=False)
    dose = np.array(dose, copy=False)

    to_minimise = create_to_minimise(net_od, dose)
    result = basinhopping(to_minimise, [np.max(dose) / np.max(net_od), 1, 1])

    return create_cal_fit(*result.x)

array([[0.  , 0.07, 0.43, 0.5 , 0.43, 0.07, 0.  ],
           [0.  , 0.14, 0.86, 1.  , 0.86, 0.14, 0.  ],
           [0.14, 0.86, 1.  , 1.  , 1.  , 0.86, 0.14],
           [0.03, 0.17, 0.2 , 0.2 , 0.2 , 0.17, 0.03]])
    """

    leaf_pair_widths = np.array(leaf_pair_widths)
    leaf_division = leaf_pair_widths / grid_resolution

    if not np.all(leaf_division.astype(int) == leaf_division):
        raise ValueError(
            "The grid resolution needs to exactly divide every leaf pair " "width."
        )

    if (
        not np.max(np.abs(jaw))  # pylint: disable = unneeded-not
        <= np.sum(leaf_pair_widths) / 2
    ):
        raise ValueError(
            "The jaw should not travel further out than the maximum leaf " "limits."
        )

    (grid, grid_leaf_map, mlc) = _determine_calc_grid_and_adjustments(
        mlc, jaw, leaf_pair_widths, grid_resolution
    )

    positions = {
        "mlc": {
            1: (-mlc[0, :, 0], -mlc[1, :, 0]),  # left
            -1: (mlc[0, :, 1], mlc[1, :, 1]),  # right
        },
        "jaw": {

def folder_analyze(volume):
    for slice in range(0, volume.shape[2]):
        stack1 = np.sum(volume[:, :, slice], axis=0)
        maxstack1 = np.max(stack1)

        stack2 = np.sum(volume[:, :, slice], axis=1)
        maxstack2 = np.max(stack2)

        if maxstack2 / maxstack1 > 1.5:  # It is a Y field folder
            return 2
        elif maxstack2 / maxstack1 < 0.5:  # It is a X field folder
            return 1
        else:
            return 3  # It is a field rotation folder

The perimeter/area data points for the relevant applicator, energy and
        ssd.
    factor_data : np.ndarray
        The insert factor data points for the relevant applicator, energy and
        ssd.

    Returns
    -------
    result : np.ndarray
        The interpolated electron insert factors for width_test and
        ratio_perim_area_test.

    """
    bbox = [
        np.min([np.min(width_data), np.min(width_test)]),
        np.max([np.max(width_data), np.max(width_test)]),
        np.min([np.min(ratio_perim_area_data), np.min(ratio_perim_area_test)]),
        np.max([np.max(ratio_perim_area_data), np.max(ratio_perim_area_test)]),
    ]

    spline = scipy.interpolate.SmoothBivariateSpline(
        width_data, ratio_perim_area_data, factor_data, kx=2, ky=1, bbox=bbox
    )

    return spline.ev(width_test, ratio_perim_area_test)

def _calc_time_steps(positions, grid_resolution, min_step_per_pixel):
    maximum_travel = []
    for _, value in positions.items():
        for _, (start, end) in value.items():
            maximum_travel.append(np.max(np.abs(end - start)))

    maximum_travel = np.max(maximum_travel)
    number_of_pixels = np.ceil(maximum_travel / grid_resolution)
    time_steps = number_of_pixels * min_step_per_pixel

    if time_steps < 10:
        time_steps = 10

    return time_steps

def _calc_time_steps(positions, grid_resolution, min_step_per_pixel):
    maximum_travel = []
    for _, value in positions.items():
        for _, (start, end) in value.items():
            maximum_travel.append(np.max(np.abs(end - start)))

    maximum_travel = np.max(maximum_travel)
    number_of_pixels = np.ceil(maximum_travel / grid_resolution)
    time_steps = number_of_pixels * min_step_per_pixel

    if time_steps < 10:
        time_steps = 10

    return time_steps

grid = dict()
    grid["jaw"] = np.arange(
        bot_grid_pos, top_grid_pos + grid_resolution, grid_resolution
    ).astype("float")

    grid_leaf_map = np.argmin(
        np.abs(grid["jaw"][:, None] - leaf_centres[None, :]), axis=1
    )

    adjusted_grid_leaf_map = grid_leaf_map - np.min(grid_leaf_map)

    leaves_to_be_calced = np.unique(grid_leaf_map)
    adjusted_mlc = mlc[:, leaves_to_be_calced, :]

    min_x = np.round(np.min(-adjusted_mlc[:, :, 0]) / grid_resolution) * grid_resolution
    max_x = np.round(np.max(adjusted_mlc[:, :, 1]) / grid_resolution) * grid_resolution

    grid["mlc"] = np.arange(min_x, max_x + grid_resolution, grid_resolution).astype(
        "float"
    )

    return grid, adjusted_grid_leaf_map, adjusted_mlc

def check_aspect_ratio(edge_lengths):
    if not np.allclose(*edge_lengths):
        if np.min(edge_lengths) > 0.95 * np.max(edge_lengths):
            raise ValueError(
                "For non-square rectangular fields, "
                "to accurately determine the rotation, "

relative_dose, depth=None, normalisation_depth=None, smoothed_normalisation=False
):
    """Normalise a pdd at a given depth. If normalisation_depth is left
    undefined then the depth of dose maximum is used for the normalisation
    depth.
    """

    if smoothed_normalisation:
        filtered = savgol_filter(relative_dose, 21, 2)
    else:
        filtered = relative_dose

    # normalisation_depth will be None if the user does not define it, if that
    # is the case simply define normalisation by 100 / the maximum value
    if normalisation_depth is None:
        normalisation = 100 / np.max(filtered)

    # However if the user did define a normalisation depth then need to
    # interpolate using the provided depth variable to find the relative dose
    # value to normalise to
    else:
        if depth is None:
            raise Exception(
                "distance variable needs to be defined to normalise to a " "depth"
            )
        interpolation = interp1d(depth, filtered)
        normalisation = 100 / interpolation(normalisation_depth)

    return relative_dose * normalisation