Я пытаюсь воспроизвести поведение следующего кода MATLAB в Python:
% Matlab code
wavelength = 10
orientation = 45
image = imread('filename.tif') % grayscale image
[mag,phase] = imgaborfilt(image, wavelength, orientation)
gabor_im = mag .* sin(phase)
К сожалению, у меня нет лицензии и я не могу запустить код. Кроме того, официальная документация Matlab для imgaborfilt не указывает точно, что делают функции.
Из-за отсутствия очевидной альтернативы я пытаюсь использовать OpenCV в Python (открыто для другие предложения). У меня нет опыта работы с OpenCV. Я пытаюсь использовать cv2.getGaborKernel и cv2.filter2D. Я также не смог найти подробную документацию о поведении этих функций. Afaik нет официальной документации обёртки Python для OpenCV. Строки документации функций предоставляют некоторую информацию, но она неполная и неточная.
Я нашел этот вопрос , где OpenCV используется в C ++ для решения проблемы. Я предполагаю, что функции работают очень похожим образом (также обратите внимание на официальную документацию C ++ ). Однако у них есть ряд дополнительных параметров. Как я могу узнать, что на самом деле функции matlab воспроизводят поведение?
# python 3.6
import numpy as np
import cv2
wavelength = 10
orientation = 45
shape = (500, 400) # arbitrary values to get running example code...
sigma = 100 # what to put for Matlab behaviour?
gamma = 1 # what to put for Matlab behaviour?
gabor_filter = cv2.getGaborKernel(shape, sigma, orientation, wavelength, gamma)
print(gabor_filter.shape) # =(401, 501). Why flipped?
image = np.random.random(shape) # create some random data.
out_buffer = np.zeros(shape)
destination_depth = -1 # like dtype for filter2D. Apparantly, -1="same as input".
thing = cv2.filter2D(image, destination_depth, gabor_filter, out_buffer)
print(out_buffer.shape, out_buffer.dtype, out_buffer.max()) # =(500, 400) float64 65.2..
print(thing.shape, thing.dtype, thing.max()) # =(500, 400) float64 65.2..
РЕДАКТИРОВАТЬ:
После получения отличного ответа Крисом Луеном go, я использовал его для создания двух функций, используя OpenCV и scikit-image соответственно, чтобы (надеюсь) воспроизвести поведение функции imgaborfit в MATLAB. Я включаю их здесь. Обратите внимание, что реализация Scikit намного медленнее, чем OpenCV.
У меня все еще есть дополнительные вопросы об этих функциях:
- Насколько точно результаты решения OpenCV и решения MATLAB согласна?
- Для людей, не желающих использовать OpenCV, я также включаю здесь решение scikit-image. Я нашел параметры, такие, что величины почти равны. Тем не менее, похоже, что фаза решения «scikit-image» отличается от решения OpenCV. Почему это?
import numpy as np
import math
import cv2
def gaborfilt_OpenCV_likeMATLAB(image, wavelength, orientation, SpatialFrequencyBandwidth=1, SpatialAspectRatio=0.5):
"""Reproduces (to what accuracy in what MATLAB version??? todo TEST THIS!) the behaviour of MATLAB imgaborfilt function using OpenCV."""
orientation = -orientation / 180 * math.pi # for OpenCV need radian, and runs in opposite direction
sigma = 0.5 * wavelength * SpatialFrequencyBandwidth
gamma = SpatialAspectRatio
shape = 1 + 2 * math.ceil(4 * sigma) # smaller cutoff is possible for speed
shape = (shape, shape)
gabor_filter_real = cv2.getGaborKernel(shape, sigma, orientation, wavelength, gamma, psi=0)
gabor_filter_imag = cv2.getGaborKernel(shape, sigma, orientation, wavelength, gamma, psi=math.pi / 2)
filtered_image = cv2.filter2D(image, -1, gabor_filter_real) + 1j * cv2.filter2D(image, -1, gabor_filter_imag)
mag = np.abs(filtered_image)
phase = np.angle(filtered_image)
return mag, phase
import numpy as np
import math
from skimage.filters import gabor
def gaborfilt_skimage_likeMATLAB(image, wavelength, orientation, SpatialFrequencyBandwidth=1, SpatialAspectRatio=0.5):
"""TODO (does not quite) reproduce the behaviour of MATLAB imgaborfilt function using skimage."""
sigma = 0.5 * wavelength * SpatialFrequencyBandwidth
filtered_image_re, filtered_image_im = gabor(
image, frequency=1 / wavelength, theta=-orientation / 180 * math.pi,
sigma_x=sigma, sigma_y=sigma/SpatialAspectRatio, n_stds=5,
)
full_image = filtered_image_re + 1j * filtered_image_im
mag = np.abs(full_image)
phase = np.angle(full_image)
return mag, phase
Код для проверки вышеуказанных функций:
from matplotlib import pyplot as plt
import numpy as np
def show(im, title=""):
plt.figure()
plt.imshow(im)
plt.title(f"{title}: dtype={im.dtype}, shape={im.shape},\n max={im.max():.3e}, min= {im.min():.3e}")
plt.colorbar()
image = np.zeros((400, 400))
image[200, 200] = 1 # a delta impulse image to visualize the filtering kernel
wavelength = 10
orientation = 33 # in degrees (for MATLAB)
mag_cv, phase_cv = gaborfilt_OpenCV_likeMATLAB(image, wavelength, orientation)
show(mag_cv, "mag") # normalized by maximum, non-zero noise even outside filter window region
show(phase_cv, "phase") # all over the place
mag_sk, phase_sk = gaborfilt_skimage_likeMATLAB(image, wavelength, orientation)
show(mag_sk, "mag skimage") # small values, zero outside filter region
show(phase_sk, "phase skimage") # and hence non-zero only inside filter window region
show(mag_cv - mag_sk/mag_sk.max(), "cv - normalized(sk)") # approximately zero-image.
show(phase_sk - phase_cv, "phase_sk - phase_cv") # phases do not agree at all! Not even in the window region!
plt.show()