python - 使用过滤器内核进行斑点检测

Question

我有高分辨率图像，其中有非常小的斑点要检测。图像的一个子部分如下所示，充满了像斑点一样的小点。 斑点的大小总是大约相同，直径为 4-5 像素。我想用像这样的手工设计的内核来捕捉它们

kernel= np.array([[ -2, -2,-1,-2,-2,-2,-1,-1,-2, -2],
                  [ -2,  0, 0, 0, 0, 0, 0, 0, 0, -2],
                  [ -2,  0, 0, 0, 0, 0, 0, 0, 0, -1],
                  [ -2,  0, 0, 2, 2, 2, 2, 0, 0, -2],
                  [ -2,  0, 0, 2, 2, 2, 2, 0, 0, -2],
                  [ -2,  0, 0, 2, 2, 2, 2, 0, 0, -2],
                  [ -2,  0, 0, 2, 2, 2, 2, 0, 0, -1],
                  [ -2,  0, 0, 0, 0, 0, 0, 0, 0, -1],
                  [ -2,  0, 0, 0, 0, 0, 0, 0, 0, -2],
                  [ -2, -2,-1,-2,-2,-2,-1,-1,-2, -2]])

并在图像上运行这个内核。

到目前为止我的代码是这样的

base_img = cv.imread(os.path.join(prepath, 'blevel.jpg'),0)

for i in range(0, base_img.shape[0]-10, 10):
    for j in range(0, base_img.shape[1]-10, 10):
        print(i,j)
        baseimg_sub = base_img[i:i+10, j:j+10]
        baseimg_sub = baseimg_sub - np.min(baseimg_sub)
        baseimg_sub = np.round(baseimg_sub / np.max(baseimg_sub), 1)
        corr_factor = np.sum(np.multiply(baseimg_sub, kernel)) / 100
        if corr_factor > 0.1 :
            cv.circle(base_img, (j+5, i+5), 7, (255,255,255))

但我抓不到他们。这些是检测到的 blob： 这似乎主要是在捕获空空间。有人可以告诉我我的内核可以改变什么吗？

标签： pythonimage-processingblobsignal-processing

Answer 1

在这种情况下，卷积似乎不可靠。我尝试了 FFT 并过滤了高频，然后简单地检测了阈值（这可以改进：注意您使用的步幅：始终是完整的内核维度；这适用于较小的步幅，例如 8；步幅较小存在风险重叠）。也许也可以使用反向转换的图像作为卷积的基础。您可能想尝试不同的过滤器值，甚至分析原始图像的光谱。这些斑点在视觉上似乎很常见。平滑的参考在代码中给出，不是我的）：

# https://akshaysin.github.io/fourier_transform.html

import numpy as np
import cv2 as cv
import os
from matplotlib import pyplot as plt
from scipy.signal import argrelextrema
from scipy import fftpack

def smooth(x, window_len=11, window='hanning'):
    """smooth the data using a window with requested size.
    
    This method is based on the convolution of a scaled window with the signal.
    The signal is prepared by introducing reflected copies of the signal 
    (with the window size) in both ends so that transient parts are minimized
    in the begining and end part of the output signal.
    
    input:
        x: the input signal 
        window_len: the dimension of the smoothing window; should be an odd integer
        window: the type of window from 'flat', 'hanning', 'hamming', 'bartlett', 'blackman'
            flat window will produce a moving average smoothing.

    output:
        the smoothed signal
        
    example:

    t=linspace(-2,2,0.1)
    x=sin(t)+randn(len(t))*0.1
    y=smooth(x)
    
    see also: 
    
    numpy.hanning, numpy.hamming, numpy.bartlett, numpy.blackman, numpy.convolve
    scipy.signal.lfilter
 
    TODO: the window parameter could be the window itself if an array instead of a string
    NOTE: length(output) != length(input), to correct this: return y[(window_len/2-1):-(window_len/2)] instead of just y.
    """
    if x.ndim != 1:
        raise ValueError # "smooth only accepts 1 dimension arrays."
    if x.size < window_len:
        raise ValueError # "Input vector needs to be bigger than window size."
    if window_len<3:
        return x
    if not window in ['flat', 'hanning', 'hamming', 'bartlett', 'blackman']:
        raise ValueError # "Window is on of 'flat', 'hanning', 'hamming', 'bartlett', 'blackman'"

    s=np.r_[x[window_len-1:0:-1],x,x[-2:-window_len-1:-1]]
    #print(len(s))
    if window == 'flat': #moving average
        w=np.ones(window_len,'d')
    else:
        w=eval('np.'+window+'(window_len)')
    y=np.convolve(w/w.sum(), s, mode='valid')
    return y

if __name__ == '__main__':
    img = cv2.imread('blevel.jpg', 0)
    dft = cv2.dft(np.float32(img), flags=cv2.DFT_COMPLEX_OUTPUT)
    dft_shift = np.fft.fftshift(dft)

    magnitude_spectrum = 200 * np.log(cv2.magnitude(dft_shift[:, :, 0], dft_shift[:, :, 1]))

    rows, cols = img.shape
    crow, ccol = int(rows / 2), int(cols / 2)  # center

    # Circular HPF mask, center circle is 0, remaining all ones
    mask = np.ones((rows, cols, 2), np.uint8)
    r = 30
    center = [crow, ccol]
    x, y = np.ogrid[:rows, :cols]
    mask_area = (x - center[0]) ** 2 + (y - center[1]) ** 2 <= r*r # circular
    #mask_area = abs((x - center[0]) + y * 0.)  <= r # directional x
    #mask_area = abs(x * 0. + (y - center[1]))  <= r # directional y
    mask[mask_area] = 0

    # apply mask and inverse DFT
    fshift = dft_shift * mask

    fshift_mask_mag = 1000 * np.log(cv2.magnitude(fshift[:, :, 0], fshift[:, :, 1]))

    f_ishift = np.fft.ifftshift(fshift)
    img_back = cv2.idft(f_ishift)
    img_back = cv2.magnitude(img_back[:, :, 0], img_back[:, :, 1])
    
    prepath = ''
    base_img = cv.imread(os.path.join(prepath, 'blevel.jpg'),0)

    plt.figure(num=None, figsize=(14, 12), dpi=80, facecolor='w', edgecolor='k')

    plt.subplot(4, 2, 1), plt.imshow(img, cmap='gray')
    plt.title('Input Image'), plt.xticks([]), plt.yticks([])
    plt.subplot(4, 2, 2), plt.imshow(magnitude_spectrum, cmap='gray')
    plt.title('After FFT'), plt.xticks([]), plt.yticks([])
    plt.subplot(4, 2, 3), plt.imshow(fshift_mask_mag, cmap='gray')
    plt.title('FFT + Mask'), plt.xticks([]), plt.yticks([])
    plt.figure(figsize=(60,40))
    plt.subplot(4, 2, 4), plt.imshow(img_back, cmap='gray')
    plt.title('After FFT Inverse'), plt.xticks([]), plt.yticks([])
    #plt.show()

    kernel= np.array([[ -2, -2,-1,-2,-2,-2,-1,-1,-2, -2],
                      [ -2,  0, 0, 0, 0, 0, 0, 0, 0, -2],
                      [ -2,  0, 0, 0, 0, 0, 0, 0, 0, -1],
                      [ -2,  0, 0, 2, 2, 2, 2, 0, 0, -2],
                      [ -2,  0, 0, 2, 2, 2, 2, 0, 0, -2],
                      [ -2,  0, 0, 2, 2, 2, 2, 0, 0, -2],
                      [ -2,  0, 0, 2, 2, 2, 2, 0, 0, -1],
                      [ -2,  0, 0, 0, 0, 0, 0, 0, 0, -1],
                      [ -2,  0, 0, 0, 0, 0, 0, 0, 0, -2],
                      [ -2, -2,-1,-2,-2,-2,-1,-1,-2, -2]])

    prepath = ''
    #base_img = cv.imread(os.path.join(prepath, 'blevel.jpg'),0)
    #blob_img = cv.imread(os.path.join(prepath, 'blob0.jpg'),0)
    #kernel = np.array(blob_img)
    #kernel = kernel - np.min(kernel)
    #kx, ky = np.shape(kernel)
    kx, ky = 8, 8 # use w/o kernel
    ksum = np.sum(kernel)
    #imgsum = np.sum(img_back)
    #print(kx, ky, kernel)
    
    invlim = np.average(img_back) + 1.5*np.std(img_back)

    for i in range(0, img_back.shape[0]-kx, int(kx/1)): # note stride
        for j in range(0, img_back.shape[1]-ky, int(ky/1)):
            #print(i,j)
            baseimg_sub = img_back[i:i+kx, j:j+ky]
            #baseimg_sub = baseimg_sub - np.min(baseimg_sub)
            #baseimg_sub = np.round(baseimg_sub / np.max(baseimg_sub), 1)
            #corr_factor = np.sum(np.multiply(baseimg_sub, kernel)) / ksum / np.sum(baseimg_sub)
            corr_factor = np.average(baseimg_sub)
            if corr_factor > invlim:
                cv.circle(base_img, (j+5, i+5), 7, (255,255,255))

    plt.figure(figsize=(60,40))
    plt.subplot(4, 2, 5), plt.imshow(base_img, cmap='gray')
    plt.show()

这会产生（注意顶部边框的一些伪影）：

这是过滤后的图像，供参考。斑点更加突出。