def add_samples_to_volume(vol, pos_ms, val):
    """
    add samples, average multiple values per voxel
    :param vol:
    :param pos_ms:
    :param val:
    :return:
    """

    # get distance between samples and their corresponding voxel centers
    pos_vs = model_space_to_volume_space(pos_ms, vol.shape[0])
    # grid_cell_centers_ms = volume_space_to_model_space(pos_vs, vol.shape[0])
    grid_cell_centers_ms = pos_ms
    dist_pos_cell_center = utils.cartesian_dist(pos_ms, grid_cell_centers_ms)

    # cluster by voxel
    unique_grid_pos, unique_counts = np.unique(pos_vs, return_counts=True, axis=0)
    values_per_voxel = np.split(val, np.cumsum(unique_counts)[:-1])
    dist_pos_cell_center_per_voxel = np.split(dist_pos_cell_center, np.cumsum(unique_counts)[:-1])
    coordinates_per_voxel = np.split(pos_vs, np.cumsum(unique_counts)[:-1])
    coordinates_per_voxel_first = np.array([c[0] for c in coordinates_per_voxel])

    # get sample closest to voxel center
    dist_pos_cell_center_per_voxel_arg_min = [np.argmin(voxel_data) for voxel_data in dist_pos_cell_center_per_voxel]
    # values_per_voxel_mean = np.array([v.mean() for v in values_per_voxel])  # TODO: weight by distance to voxel center
    values_per_voxel_closest = np.array([v[dist_pos_cell_center_per_voxel_arg_min[vi]]
                                        for vi, v in enumerate(values_per_voxel)])
    vol[coordinates_per_voxel_first[:, 0], coordinates_per_voxel_first[:, 1], coordinates_per_voxel_first[:, 2]] = \
        values_per_voxel_closest  # values_per_voxel_mean
    return vol

    add samples, average multiple values per voxel
    :param vol:
    :param pos_ms:
    :param val:
    :return:
    
source/sdf.py