Python skimage.feature.blob_dog() Examples

def __init__(
            self,
            min_sigma: Union[Number, Tuple[Number, ...]],
            max_sigma: Union[Number, Tuple[Number, ...]],
            num_sigma: int,
            threshold: Number,
            overlap: float = 0.5,
            measurement_type='max',
            is_volume: bool = True,
            detector_method: str = 'blob_log',
            exclude_border: Optional[int] = None,
    ) -> None:

        self.min_sigma = min_sigma
        self.max_sigma = max_sigma
        self.num_sigma = num_sigma
        self.threshold = threshold
        self.overlap = overlap
        self.is_volume = is_volume
        self.measurement_function = self._get_measurement_function(measurement_type)
        self.exclude_border = exclude_border
        try:
            self.detector_method = blob_detectors[detector_method]
        except ValueError:
            raise ValueError("Detector method must be one of {blob_log, blob_dog, blob_doh}") 

def test_dog_feature(self):
        blob = self.patch.data['blob_dog']
        delta = 1e-4

        test_min = np.min(blob)
        exp_min = 0.0
        self.assertAlmostEqual(test_min, exp_min, delta=delta, msg="Expected min {}, got {}".format(exp_min, test_min))

        test_max = np.max(blob)
        exp_max = 37.9625
        self.assertAlmostEqual(test_max, exp_max, delta=delta, msg="Expected max {}, got {}".format(exp_max, test_max))

        test_mean = np.mean(blob)
        exp_mean = 0.0545
        self.assertAlmostEqual(test_mean, exp_mean, delta=delta,
                               msg="Expected mean {}, got {}".format(exp_mean, test_mean))

        test_median = np.median(blob)
        exp_median = 0.0
        self.assertAlmostEqual(test_median, exp_median, delta=delta,
                               msg="Expected median {}, got {}".format(exp_median, test_median)) 

def _peak_find_dog_chunk(data, **kwargs):
    """Find peaks in a chunk using skimage's blob_dog function.

    Parameters
    ----------
    data : NumPy array
    min_sigma : float, optional
    max_sigma : float, optional
    sigma_ratio : float, optional
    threshold : float, optional
    overlap : float, optional
    normalize_value : float, optional
        All the values in data will be divided by this value.
        If no value is specified, the max value in each individual image will
        be used.

    Returns
    -------
    peak_array : NumPy 2D object array
        Same size as the two last dimensions in data.
        The peak positions themselves are stored in 2D NumPy arrays
        inside each position in peak_array. This is done instead of
        making a 4D NumPy array, since the number of found peaks can
        vary in each position.

    Example
    -------
    >>> s = pxm.dummy_data.dummy_data.get_cbed_signal()
    >>> import pyxem.utils.dask_tools as dt
    >>> peak_array = dt._peak_find_dog_chunk(s.data)
    >>> peaks00 = peak_array[0, 0]
    >>> peaks23 = peak_array[2, 3]

    """
    output_array = np.empty(data.shape[:-2], dtype="object")
    for index in np.ndindex(data.shape[:-2]):
        islice = np.s_[index]
        output_array[islice] = _peak_find_dog_single_frame(image=data[islice], **kwargs)
    return output_array 

def setUpClass(cls):
        cls.patch = EOPatch.load(cls.TEST_PATCH_FILENAME)
        cls._prepare_patch(cls.patch)

        BlobTask((FeatureType.DATA, 'ndvi', 'blob'), blob_dog, sigma_ratio=1.6, min_sigma=1, max_sigma=30,
                 overlap=0.5, threshold=0).execute(cls.patch)
        DoGBlobTask((FeatureType.DATA, 'ndvi', 'blob_dog'), threshold=0).execute(cls.patch)
        LoGBlobTask((FeatureType.DATA, 'ndvi', 'blob_log'), log_scale=True, threshold=0).execute(cls.patch)
        DoHBlobTask((FeatureType.DATA, 'ndvi', 'blob_doh'), num_sigma=5, threshold=0).execute(cls.patch)

        cls.initial_patch = EOPatch.load(cls.TEST_PATCH_FILENAME)
        cls._prepare_patch(cls.initial_patch) 

def test_blob_feature(self):
        self.assertTrue(np.allclose(self.patch.data['blob'], self.patch.data['blob_dog']),
                        msg='DoG derived class result not equal to base class result') 

def blob_image_multiscale2(image, type=0,scale=2):
    # function that return a list of blob_coordinates, 0 = dog, 1 = doh, 2 = log
    list = []
    image = norm.normalize(image)
    for z, slice in tqdm(enumerate(image)):
        # init list of different sigma/zoom blobs
        featureblobs = []
        # x = 0,1,2,3,4
        if scale == 2:
            # for x in xrange(0,6):
            #     if type == 0:
            #         featureblobs.append(feature.blob_dog(slice, 2**x, 2**x))
            #     if type == 1:
            #         featureblobs.append(feature.blob_doh(slice, 2**x, 2**x))
            #     if type == 2:
            #         featureblobs.append(feature.blob_log(slice, 2**x, 2**x))
            for x in xrange(0,5):
                if type == 0:
                    featureblobs.append(feature.blob_dog(slice, 2**x, 2**(x+1)))
                if type == 1:
                    featureblobs.append(feature.blob_doh(slice, 2**x, 2**(x+1)))
                if type == 2:
                    featureblobs.append(feature.blob_log(slice, 2**x, 2**(x+1),16,.1))
        else:
            for x in xrange(0,4):
                if type == 0:
                    featureblobs.append(feature.blob_dog(slice, 3**x, 3**x))
                if type == 1:
                    featureblobs.append(feature.blob_doh(slice, 3**x, 3**x))
                if type == 2:
                    featureblobs.append(feature.blob_log(slice, 3**x, 3**x))
        # init list of blob coords
        blob_coords = []
        #print featureblobs
        # start at biggest blob size
        for featureblob in reversed(featureblobs):
            # for every blob found of a blobsize

            for blob in enumerate(featureblob):
                # if that blob is not within range of another blob, add it
                blob = blob[1]
                if not within_range(blob, blob_coords):
                    blob_coords.append([z, blob[0], blob[1], blob[2]])
        list.append(blob_coords)
    return list 

def find_peaks_dog(
    z,
    min_sigma=1.0,
    max_sigma=50.0,
    sigma_ratio=1.6,
    threshold=0.2,
    overlap=0.5,
    exclude_border=False,
):
    """
    Finds peaks via the difference of Gaussian Matrices method from
    `scikit-image`.

    Parameters
    ----------
    z : numpy.ndarray
        2-d array of intensities
    float min_sigma, max_sigma, sigma_ratio, threshold, overlap
        Additional parameters to be passed to the algorithm. See `blob_dog`
        documentation for details:
        http://scikit-image.org/docs/dev/api/skimage.feature.html#blob-dog

    Returns
    -------
    peaks : numpy.ndarray
        Array of peak pixel coordinates with shape (n_peaks, 2).

    Notes
    -----
    While highly effective at finding even very faint peaks, this method is
    sensitive to fluctuations in intensity near the edges of the image.

    """
    z = z / np.max(z)
    blobs = blob_dog(
        z,
        min_sigma=min_sigma,
        max_sigma=max_sigma,
        sigma_ratio=sigma_ratio,
        threshold=threshold,
        overlap=overlap,
    )

    centers = blobs[:, :2]
    return centers 

def _peak_find_dog_single_frame(
    image,
    min_sigma=0.98,
    max_sigma=55,
    sigma_ratio=1.76,
    threshold=0.36,
    overlap=0.81,
    normalize_value=None,
):
    """Find peaks in a single frame using skimage's blob_dog function.

    Parameters
    ----------
    image : NumPy 2D array
    min_sigma : float, optional
    max_sigma : float, optional
    sigma_ratio : float, optional
    threshold : float, optional
    overlap : float, optional
    normalize_value : float, optional
        All the values in image will be divided by this value.
        If no value is specified, the max value in the image will be used.

    Returns
    -------
    peaks : NumPy 2D array
        In the form [[x0, y0], [x1, y1], [x2, y2], ...]

    Example
    -------
    >>> s = pxm.dummy_data.dummy_data.get_cbed_signal()
    >>> import pyxem.utils.dask_tools as dt
    >>> peaks = _peak_find_dog_single_frame(s.data[0, 0])

    """

    if normalize_value is None:
        normalize_value = np.max(image)
    peaks = blob_dog(
        image / normalize_value,
        min_sigma=min_sigma,
        max_sigma=max_sigma,
        sigma_ratio=sigma_ratio,
        threshold=threshold,
        overlap=overlap,
        exclude_border=False,
    )
    peak = peaks[:, :2].astype(np.float64)

    return peak 

def _peak_find_dog(dask_array, **kwargs):
    """Find peaks in a dask array using skimage's blob_dog function.

    Parameters
    ----------
    dask_array : Dask array
        Must be at least 2 dimensions.
    min_sigma : float, optional
    max_sigma : float, optional
    sigma_ratio : float, optional
    threshold : float, optional
    overlap : float, optional
    normalize_value : float, optional
        All the values in dask_array will be divided by this value.
        If no value is specified, the max value in each individual image will
        be used.

    Returns
    -------
    peak_array : dask object array
        Same size as the two last dimensions in data.
        The peak positions themselves are stored in 2D NumPy arrays
        inside each position in peak_array. This is done instead of
        making a 4D NumPy array, since the number of found peaks can
        vary in each position.

    Example
    -------
    >>> s = pxm.dummy_data.dummy_data.get_cbed_signal()
    >>> import dask.array as da
    >>> dask_array = da.from_array(s.data, chunks=(5, 5, 25, 25))
    >>> import pyxem.utils.dask_tools as dt
    >>> peak_array = _peak_find_dog(dask_array)
    >>> peak_array_computed = peak_array.compute()

    """
    dask_array_rechunked = _rechunk_signal2d_dim_one_chunk(dask_array)
    drop_axis = (dask_array_rechunked.ndim - 2, dask_array_rechunked.ndim - 1)
    output_array = da.map_blocks(
        _peak_find_dog_chunk,
        dask_array_rechunked,
        drop_axis=drop_axis,
        dtype=np.object,
        **kwargs
    )
    return output_array 

def blob_detection(
    input_image,
    method="log",
    threshold=0.5,
    min_sigma=3,
    max_sigma=20,
    overlap=0.5,
    return_sigmas=False,
    **kwargs,
):
    """
    Interface to the `feature.blob_*`_ methods implemented in scikit-image. A
    blob is defined as local a maximum of a Gaussian-filtered image.

    .. _`feature.blob_*`:\
    https://scikit-image.org/docs/dev/auto_examples/features_detection/plot_blob.html

    Parameters
    ----------
    input_image : array_like
        Array of shape (m, n) containing the input image. Nan values are ignored.
    method : {'log', 'dog', 'doh'}, optional
        The method to use: 'log' = Laplacian of Gaussian, 'dog' = Difference of
        Gaussian, ''
    threshold : float, optional
        Detection threshold.
    min_sigma : float, optional
        The minimum standard deviation for the Gaussian kernel.
    max_sigma : float, optional
        The maximum standard deviation for the Gaussian kernel.
    overlap : float, optional
        A value between 0 and 1. If the area of two blobs overlaps by a fraction
        greater than threshold, the smaller blob is eliminated.
    return_sigmas : bool, optional
        If True, the return array has a third column indicating the standard
        deviations of the Gaussian kernels that detected the blobs.
    """
    if method not in ["log", "dog", "doh"]:
        raise ValueError("unknown method %s, must be 'log', 'dog' or 'doh'" % method)

    if method == "log":
        detector = feature.blob_log
    elif method == "dog":
        detector = feature.blob_dog
    else:
        detector = feature.blob_doh

    blobs = detector(
        input_image,
        min_sigma=min_sigma,
        max_sigma=max_sigma,
        threshold=threshold,
        overlap=overlap,
        **kwargs,
    )

    if not return_sigmas:
        blobs = blobs[:, :2]

    return np.column_stack([blobs[:, 1], blobs[:, 0]])