Harris Corner Detector生成的特征对于缩放不是不变的。对于特征跟踪,我们需要对仿射变换不变的特征。拉普拉斯斑点检测器是生成对缩放不变的特征的基本方法之一。
拉普拉斯斑点检测器的想法是将图像与多个尺度的“斑点滤波器”进行卷积,并在得到的尺度空间中寻找滤波器响应的极值。
斑点滤波器(Blob Filter):
这个滤波器是由沿x轴和y轴的double derivating Gaussian filter产生并相加。以上也称为高斯拉普拉斯算子。
当图像中有一个角时,LOG滤波器的输出最大。
我们想通过将图像与不同的sigma卷积来找到blob的特征尺度。但随着我们增加sigma,blob滤波器对图像的响应会降低,因此我们可能无法找到更好的尺度(scale)。为了减少这种影响,我们将LOG滤波器σ平方相乘。这个过程称为尺度归一化。
拉普拉斯档案对半径为r的binary圆的最大响应是σ= 1.414 * r
以上是blob滤波器的一些基础知识。整个过程归结为两个步骤
- 在几个尺度上用尺度归一化拉普拉斯变换卷积图像(不同尺度意味着不同的西格玛)
- 在尺度空间中找到平方拉普拉斯响应的最大值
拉普拉斯Blob检测的Python代码
该代码有四个主要步骤:
- 生成LOG滤波器
- 用高斯滤波器对图像进行卷积
- 找到最大峰值
- 绘制blob。
生成LOG滤波器的Python代码
import cv2 from pylab import * import numpy as np import matplotlib.pyplot as plt img = cv2.imread("sunflowers.jpg",0) #gray scale conversion from scipy import ndimage from scipy.ndimage import filters from scipy import spatial k = 1.414 sigma = 1.0 img = img/255.0 #image normalization def LoG(sigma): #window size n = np.ceil(sigma*6) y,x = np.ogrid[-n//2:n//2+1,-n//2:n//2+1] y_filter = np.exp(-(y*y/(2.*sigma*sigma))) x_filter = np.exp(-(x*x/(2.*sigma*sigma))) final_filter = (-(2*sigma**2) + (x*x + y*y) ) * (x_filter*y_filter) * (1/(2*np.pi*sigma**4)) return final_filter
LoG函数将sigma作为输入,并生成大小为6 * sigma(最佳选项)的滤波器。
用高斯滤波器对图像进行卷积
图像与不同的LOG滤波器进行卷积,Python代码如下:
def LoG_convolve(img): log_images = [] #to store responses for i in range(0,9): y = np.power(k,i) sigma_1 = sigma*y #sigma filter_log = LoG(sigma_1) #filter generation image = cv2.filter2D(img,-1,filter_log) # convolving image image = np.pad(image,((1,1),(1,1)),'constant') #padding image = np.square(image) # squaring the response log_images.append(image) log_image_np = np.array([i for i in log_images]) # storing the #in numpy array return log_image_np log_image_np = LoG_convolve(img) #print(log_image_np.shape)
要选择不同的sigma,我将sigma与k ^ i相乘。log_image_np数组的维度为z * image_height * image_width,其中z表示与sigma相乘的第k次方。
在上面的函数中,我只生成了9个滤波器,图中的9个滤波器是
找到最大峰值
def detect_blob(log_image_np): co_ordinates = [] #to store co ordinates (h,w) = img.shape for i in range(1,h): for j in range(1,w): slice_img = log_image_np[:,i-1:i+2,j-1:j+2] #9*3*3 slice result = np.amax(slice_img) #finding maximum if result >= 0.03: #threshold z,x,y = np.unravel_index(slice_img.argmax(),slice_img.shape) co_ordinates.append((i+x-1,j+y-1,k**z*sigma)) #finding co-rdinates return co_ordinates co_ordinates = list(set(detect_blob(log_image_np)))
在每个像素处,在所有尺度上考虑3x3邻域,其给出9x3x3矩阵。矩阵的最大峰值由矩阵的最大值给出。
对于每个像素,都会有一个最大值。但并非所有像素都有助于blob。所以我们需要通过阈值来消除它们。
如果最大值大于阈值,则考虑该点并存储坐标。
(z,x,y)是切片图像的坐标,因此我们需要将它们更改为原始图像坐标。在更改它们之后,我们将它们存储在co_ordinates列表中并绘制它们。Python代码如下:
fig, ax = plt.subplots() nh,nw = img.shape count = 0 ax.imshow(img, interpolation='nearest',cmap="gray") for blob in co_ordinates: y,x,r = blob c = plt.Circle((x, y), r*1.414, color='red', linewidth=1.5, fill=False) ax.add_patch(c) ax.plot() plt.show()
上面的代码在存储列表的每个点绘制blob。
结果
在上面的图像中,我们可以发现许多blob都是重叠的。为了减少重叠的斑点,我们将找到斑点之间的重叠区域。如果该区域大于阈值,则将丢弃具有较小半径的斑点。python实现如下:
#provided in scipy doucumentaion def blob_overlap(blob1, blob2): n_dim = len(blob1) - 1 root_ndim = sqrt(n_dim) #print(n_dim) # radius of two blobs r1 = blob1[-1] * root_ndim r2 = blob2[-1] * root_ndim d = sqrt(np.sum((blob1[:-1] - blob2[:-1])**2)) #no overlap between two blobs if d > r1 + r2: return 0 # one blob is inside the other, the smaller blob must die elif d <= abs(r1 - r2): return 1 else: #computing the area of overlap between blobs ratio1 = (d ** 2 + r1 ** 2 - r2 ** 2) / (2 * d * r1) ratio1 = np.clip(ratio1, -1, 1) acos1 = math.acos(ratio1) ratio2 = (d ** 2 + r2 ** 2 - r1 ** 2) / (2 * d * r2) ratio2 = np.clip(ratio2, -1, 1) acos2 = math.acos(ratio2) a = -d + r2 + r1 b = d - r2 + r1 c = d + r2 - r1 d = d + r2 + r1 area = (r1 ** 2 * acos1 + r2 ** 2 * acos2 -0.5 * sqrt(abs(a * b * c * d))) return area/(math.pi * (min(r1, r2) ** 2)) def redundancy(blobs_array, overlap): sigma = blobs_array[:, -1].max() distance = 2 * sigma * sqrt(blobs_array.shape[1] - 1) tree = spatial.cKDTree(blobs_array[:, :-1]) pairs = np.array(list(tree.query_pairs(distance))) if len(pairs) == 0: return blobs_array else: for (i, j) in pairs: blob1, blob2 = blobs_array[i], blobs_array[j] if blob_overlap(blob1, blob2) > overlap: if blob1[-1] > blob2[-1]: blob2[-1] = 0 else: blob1[-1] = 0 return np.array([b for b in blobs_array if b[-1] > 0]) co_ordinates = np.array(co_ordinates) co_ordinates = redundancy(co_ordinates,0.5)
删除重叠blob后输出为
我们可以看到大多数斑点已被消除。