python类peak_local_max()的实例源码

def circles_per_radius(hough_radii, hough_res, number_circles_per_radius=16):
    """
    Rates found circles by intensity of peaks inside the Hough matrix. 
    Chooses the best circles for each radius
    Returns: Selected circle centers with their quality and radius      
    """
    centers = []
    accums = []
    radii = []

    # for each radius and hough peak image
    # zip([32,33,34,35,36],[(32, hough_peaks_img),(33, hough_peaks_img),(34, hough_peaks_img), ... ])
    for radius, h in zip(hough_radii, hough_res): # iterieren durch circles (h)

        # sort peaks, which represent the quality of circles by intensity
        peaks = peak_local_max(h, num_peaks=number_circles_per_radius)
        centers.extend(peaks)
        accums.extend(h[peaks[:, 0], peaks[:, 1]])
        # iterate through every (y,x) in peaks and get corresponding color value from h, which represents quality value of circle (?)
        # so acuums represents quality
        radii.extend([radius] * number_circles_per_radius)
        #

    return (centers, accums, radii)

def get_echos(dat_arr,min_distance=10):
    """
    getting local maxima and minima from 1D data -> e.g. speckle echos
    strategy: using peak_local_max (from skimage) with min_distance parameter to find well defined local maxima
    using np.argmin to find absolute minima between relative maxima
    returns [max_ind,min_ind] -> lists of indices corresponding to local maxima/minima
    by LW 10/23/2018
    """
    from skimage.feature import peak_local_max
    max_ind=peak_local_max(dat_arr, min_distance)   # !!! careful, skimage function reverses the order (wtf?)
    min_ind=[]
    for i in range(len(max_ind[:-1])):
        min_ind.append(max_ind[i+1][0]+np.argmin(dat_arr[max_ind[i+1][0]:max_ind[i][0]]))
    #unfortunately, skimage function fu$$s up the format: max_ind is an array of a list of lists...fix this:
    mmax_ind=[]
    for l in max_ind:
        mmax_ind.append(l[0])
    #return [mmax_ind,min_ind]
    return [list(reversed(mmax_ind)),list(reversed(min_ind))]

def findMaximaOnFG(self, param):
        self.defineFG(param)
        #self.smooth_corr()
        self.coorExtract = [0, 0]
        xmin, ymin = self.coorExtract


        img=self.candi
        img [self.FG ==0] =0

        im = img_as_float(img)
        image_max = ndimage.maximum_filter(im, size=10, mode='constant')
        coordinates = peak_local_max(im, min_distance=10)

        tep=np.zeros(self.candi.shape)
        for i,ide in enumerate(coordinates):
            tep[ide[0],ide[1]] = self.candi[ide[0],ide[1]]
        lbl = ndimage.label(tep)[0]
        centerc = np.round(ndimage.measurements.center_of_mass(tep, lbl, range(1,np.max(lbl)+1)))
        if centerc.size > 0:
            self.centersX = centerc[:,0].astype(int)
            self.centersY = centerc[:,1].astype(int)
        self.nComponents = len(self.centersX)

def find_candidates(greyscale_image, neighborhood_size=20, candidate_threshold=.5):
    corner_likelihood = compute_inital_corner_likelihood(greyscale_image)
    # TODO: the absolute threshold should be statistically determined based on actual checkerboard images
    candidates = peak_local_max(corner_likelihood, neighborhood_size, corner_likelihood.max() * candidate_threshold)
    return candidates

def process(self, im):
        (width, height, _) = im.image.shape

        img_adapted = im.prep(self.transform)

        if width > self.max_resized or height > self.max_resized:
            scale_height = self.max_resized / height
            scale_width = self.max_resized / width
            scale = min(scale_height, scale_width)
            img_adapted = resize(img_adapted, (int(width * scale), int(height * scale)))

        edges = canny(img_adapted, sigma=self.sigma)

        # Detect two radii
        # Calculate image diameter
        shape = im.image.shape
        diam = math.sqrt(shape[0] ** 2 + shape[1] ** 2)
        radii = np.arange(diam / 3, diam * 0.8, 2)
        hough_res = hough_circle(edges, radii)

        accums = []
        for radius, h in zip(radii, hough_res):
            # For each radius, extract two circles
            peaks = peak_local_max(h, num_peaks=1, min_distance=1)
            if len(peaks) > 0:
                accums.extend(h[peaks[:, 0], peaks[:, 1]])

        if len(accums) == 0:  # TODO: fix, should not happen
            return [0]

        idx = np.argmax(accums)
        return [accums[idx]]

def makeSegmentation(d, s, nbins = 256, verbose = True, sigma = 0.75, min_distance = 1):
  fn = os.path.join(datadir, 'data_stage%d_%s.npy' % (s,features[0]));
  dists = np.load(fn);
  fn = os.path.join(datadir, 'data_stage%d_%s.npy' % (s,features[1]));
  rots = np.load(fn);

  ddata = np.vstack([np.log(dists[:,d]), (rots[:,d])]).T
  #gmmdata = np.vstack([dists[:,j], (rots[:,j])]).T
  #ddata.shape

  nanids = np.logical_or(np.any(np.isnan(ddata), axis=1), np.any(np.isinf(ddata), axis=1));
  ddata = ddata[~nanids,:];
  #ddata.shape
  imgbin = None;
  img2 = smooth(ddata, nbins = [nbins, nbins], sigma = (sigma,sigma))
  #img = smooth(ddata, nbins = [nbins, nbins], sigma = (1,1))

  local_maxi = peak_local_max(img2, indices=False, min_distance = min_distance)
  imgm2 = img2.copy();
  imgm2[local_maxi] = 3 * imgm2.max();

  if verbose:
    imgbin = smooth(ddata, nbins = [nbins, nbins], sigma = None)
    plt.figure(220); plt.clf()
    plt.subplot(2,2,1)
    plt.imshow(imgbin)
    plt.subplot(2,2,2)
    plt.imshow(img2)
    plt.subplot(2,2,3);
    plt.imshow(imgm2,  cmap=plt.cm.jet, interpolation='nearest')

  markers = ndi.label(local_maxi)[0]
  labels = watershed(-img2, markers, mask = None);
  print "max labels: %d" % labels.max()

  if verbose:
    fig, axes = plt.subplots(ncols=3, sharex=True, sharey=True, subplot_kw={'adjustable':'box-forced'})
    ax0, ax1, ax2 = axes
    ax0.imshow(img2, cmap=plt.cm.jet, interpolation='nearest')
    ax0.set_title('PDF')
    labels[imgbin==0] = 0;
    labels[0,0] = -1;
    ax1.imshow(labels, cmap=plt.cm.rainbow, interpolation='nearest')
    ax1.set_title('Segmentation on Data')
    #labelsws[0,0] = -1;
    #ax2.imshow(labelsws, cmap=plt.cm.rainbow, interpolation='nearest')
    #ax1.set_title('Segmentation Full')

  return labels;


#classification for a specific work based on the segmentation above