The convolution seems to be unreliable in this case. I have tried a FFT and filtered the high frequencies and then simply detected threshold values (this can be improved: note the stride you use: always a full kernel dimension; this works for smaller strides, e.g. 8; with smaller strides there is a risk of overlap). Maybe also use the back-transformed image as a basis for the convolution. You may want to try different filter values or even analyse the spectrum of the original image. The blobs seem to be quite frequent visually. The reference for the smooth is given in the code, not mine):
# 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()
This produces (note some artifacts at the top border):
This is the filtered image, for reference. The blobs are more prominent.
The blobs are always about the same size 4-5px diameter.
I want to catch them using a hand-engineered kernel such as this
Which seems to be catching emtpy spaces mostly. Can someone tell me please what can I change in my kernel?
