<SYSTEM_TASK:>
Solve the following problem using Python, implementing the functions described below, one line at a time
<END_TASK>
<USER_TASK:>
Description:
 def porosimetry(im, sizes=25, inlets=None, access_limited=True, mode='hybrid'):
 r""" Performs a porosimetry simulution on the image Parameters im : ND-array An ND image of the porous material containing True values in the pore space. sizes : array_like or scalar The sizes to invade. If a list of values of provided they are used directly. If a scalar is provided then that number of points spanning the min and max of the distance transform are used. inlets : ND-array, boolean A boolean mask with True values indicating where the invasion enters the image. By default all faces are considered inlets, akin to a mercury porosimetry experiment. Users can also apply solid boundaries to their image externally before passing it in, allowing for complex inlets like circular openings, etc. This argument is only used if ``access_limited`` is ``True``. access_limited : Boolean This flag indicates if the intrusion should only occur from the surfaces (``access_limited`` is True, which is the default), or if the invading phase should be allowed to appear in the core of the image. The former simulates experimental tools like mercury intrusion porosimetry, while the latter is useful for comparison to gauge the extent of shielding effects in the sample. mode : string Controls with method is used to compute the result. Options are: 'hybrid' - (default) Performs a distance tranform of the void space, thresholds to find voxels larger than ``sizes[i]``, trims the resulting mask if ``access_limitations`` is ``True``, then dilates it using the efficient fft-method to obtain the non-wetting fluid configuration. 'dt' - Same as 'hybrid', except uses a second distance transform, relative to the thresholded mask, to find the invading fluid configuration. The choice of 'dt' or 'hybrid' depends on speed, which is system and installation specific. 'mio' - Using a single morphological image opening step to obtain the invading fluid confirguration directly, *then* trims if ``access_limitations`` is ``True``. This method is not ideal and is included mostly for comparison purposes. The morphological operations are done using fft-based method implementations. Returns ------- image : ND-array A copy of ``im`` with voxel values indicating the sphere radius at which it becomes accessible from the inlets. This image can be used to find invading fluid configurations as a function of applied capillary pressure by applying a boolean comparison: ``inv_phase = im > r`` where ``r`` is the radius (in voxels) of the invading sphere. Of course, ``r`` can be converted to capillary pressure using your favorite model. Notes ----- There are many ways to perform this filter, and PoreSpy offer 3, which users can choose between via the ``mode`` argument. These methods all work in a similar way by finding which foreground voxels can accomodate a sphere of a given radius, then repeating for smaller radii. See Also -------- fftmorphology local_thickness """

    if im.ndim != im.squeeze().ndim:
        warnings.warn('Input image conains a singleton axis:' + str(im.shape) +
                      ' Reduce dimensionality with np.squeeze(im) to avoid' +
                      ' unexpected behavior.')

    dt = spim.distance_transform_edt(im > 0)

    if inlets is None:
        inlets = get_border(im.shape, mode='faces')

    if isinstance(sizes, int):
        sizes = sp.logspace(start=sp.log10(sp.amax(dt)), stop=0, num=sizes)
    else:
        sizes = sp.unique(sizes)[-1::-1]

    if im.ndim == 2:
        strel = ps_disk
    else:
        strel = ps_ball

    if mode == 'mio':
        pw = int(sp.floor(dt.max()))
        impad = sp.pad(im, mode='symmetric', pad_width=pw)
        inletspad = sp.pad(inlets, mode='symmetric', pad_width=pw)
        inlets = sp.where(inletspad)
#        sizes = sp.unique(sp.around(sizes, decimals=0).astype(int))[-1::-1]
        imresults = sp.zeros(sp.shape(impad))
        for r in tqdm(sizes):
            imtemp = fftmorphology(impad, strel(r), mode='erosion')
            if access_limited:
                imtemp = trim_disconnected_blobs(imtemp, inlets)
            imtemp = fftmorphology(imtemp, strel(r), mode='dilation')
            if sp.any(imtemp):
                imresults[(imresults == 0)*imtemp] = r
        imresults = extract_subsection(imresults, shape=im.shape)
    elif mode == 'dt':
        inlets = sp.where(inlets)
        imresults = sp.zeros(sp.shape(im))
        for r in tqdm(sizes):
            imtemp = dt >= r
            if access_limited:
                imtemp = trim_disconnected_blobs(imtemp, inlets)
            if sp.any(imtemp):
                imtemp = spim.distance_transform_edt(~imtemp) < r
                imresults[(imresults == 0)*imtemp] = r
    elif mode == 'hybrid':
        inlets = sp.where(inlets)
        imresults = sp.zeros(sp.shape(im))
        for r in tqdm(sizes):
            imtemp = dt >= r
            if access_limited:
                imtemp = trim_disconnected_blobs(imtemp, inlets)
            if sp.any(imtemp):
                imtemp = fftconvolve(imtemp, strel(r), mode='same') > 0.0001
                imresults[(imresults == 0)*imtemp] = r
    else:
        raise Exception('Unreckognized mode ' + mode)
    return imresults