实现日志 Gabor 滤波器组

问题描述 投票:0回答:2

我正在阅读这篇论文“自可逆 2D Log-Gabor 小波”它定义了 2D Log-Gabor 滤波器:

论文还指出滤波器仅覆盖频率空间的一侧,并在该图像中显示

在我尝试实现过滤器时,我得到的结果与论文中所说的不符。让我从我的实现开始,然后我将说明问题。

实施:

  1. 我创建了一个包含滤波器的二维数组,并转换了每个索引,以便频域的原点位于数组的中心,正 x 轴向右,正 y 轴向上。

    number_scales = 5         # scale resolution
    number_orientations = 9   # orientation resolution
    N = constantDim           # image dimensions
    
    def getLogGaborKernal(scale, angle, logfun=math.log2, norm = True):
        # setup up filter configuration
        center_scale = logfun(N) - scale          
        center_angle = ((np.pi/number_orientations) * angle) if (scale % 2) \
                    else ((np.pi/number_orientations) * (angle+0.5))
        scale_bandwidth =  0.996 * math.sqrt(2/3)
        angle_bandwidth =  0.996 * (1/math.sqrt(2)) * (np.pi/number_orientations)
    
        # 2d array that will hold the filter
        kernel = np.zeros((N, N))
        # get the center of the 2d array so we can shift origin
        middle = math.ceil((N/2)+0.1)-1
    
        # calculate the filter
        for x in range(0,constantDim):
            for y in range(0,constantDim):
                # get the transformed x and y where origin is at center
                # and positive x-axis goes right while positive y-axis goes up
                x_t, y_t = (x-middle),-(y-middle)
                # calculate the filter value at given index
                kernel[y,x] = logGaborValue(x_t,y_t,center_scale,center_angle,
            scale_bandwidth, angle_bandwidth,logfun)
    
        # normalize the filter energy
        if norm:
            Kernel = kernel / np.sum(kernel**2)
        return kernel
    
  2. 为了计算每个索引处的滤波器值,我们在对数极空间中进行另一个变换

    def logGaborValue(x,y,center_scale,center_angle,scale_bandwidth,
                  angle_bandwidth, logfun):
        # transform to polar coordinates
        raw, theta = getPolar(x,y)
        # if we are at the center, return 0 as in the log space
        # zero is not defined
        if raw == 0:
            return 0
    
        # go to log polar coordinates
        raw = logfun(raw)
    
        # calculate (theta-center_theta), we calculate cos(theta-center_theta) 
        # and sin(theta-center_theta) then use atan to get the required value,
        # this way we can eliminate the angular distance wrap around problem
        costheta, sintheta = math.cos(theta), math.sin(theta)
        ds = sintheta * math.cos(center_angle) - costheta * math.sin(center_angle)    
        dc = costheta * math.cos(center_angle) + sintheta * math.sin(center_angle)  
        dtheta = math.atan2(ds,dc)
    
        # final value, multiply the radial component by the angular one
        return math.exp(-0.5 * ((raw-center_scale) / scale_bandwidth)**2) * \
                math.exp(-0.5 * (dtheta/angle_bandwidth)**2)
    

问题:

  1. 角度:论文指出,对从 1->8 的角度进行索引可以很好地覆盖方向,但在我的实现中,从 1->n 的角度不会覆盖除半方向之外的情况。甚至垂直方向也没有被正确覆盖。这可以如图所示,其中包含尺度为 3 且方向范围为 1->8 的滤波器组:

  2. 覆盖范围:从上面的过滤器可以明显看出过滤器覆盖了空间的两侧,这不是论文所说的。通过使用范围从 -4 -> 4 的 9 个方向,可以使这一点更加明确。下图包含一张图像中的所有滤波器,以显示它如何覆盖光谱的两侧(该图像是通过在每个位置取最大值来创建的)来自所有过滤器):

  3. 中间列(方向 $\pi / 2$): 在方向从 3 -> 8 的第一个图中可以看出,滤波器在方向 $\pi / 2$ 处消失。这是正常的吗?当我将所有滤镜(所有 5 个尺度和 9 个方向)组合到一张图像中时,也可以看到这一点:

更新: 在空间域中添加滤波器的脉冲响应,可以看到 -4 和 4 方向有明显的失真:

python image-processing computer-vision
2个回答
6
投票

经过大量代码分析,我发现我的实现是正确的,但是

getPolar
函数被搞乱了,所以上面的代码应该可以正常工作。如果有人正在寻找的话,这是一个没有
getPolar
功能的新代码:

number_scales = 5          # scale resolution
number_orientations = 8    # orientation resolution
N = 128                    # image dimensions
def getFilter(f_0, theta_0):
    # filter configuration
    scale_bandwidth =  0.996 * math.sqrt(2/3)
    angle_bandwidth =  0.996 * (1/math.sqrt(2)) * (np.pi/number_orientations)

    # x,y grid
    extent = np.arange(-N/2, N/2 + N%2)
    x, y = np.meshgrid(extent,extent)

    mid = int(N/2)
    ## orientation component ##
    theta = np.arctan2(y,x)
    center_angle = ((np.pi/number_orientations) * theta_0) if (f_0 % 2) \
                else ((np.pi/number_orientations) * (theta_0+0.5))

    # calculate (theta-center_theta), we calculate cos(theta-center_theta) 
    # and sin(theta-center_theta) then use atan to get the required value,
    # this way we can eliminate the angular distance wrap around problem
    costheta = np.cos(theta)
    sintheta = np.sin(theta)
    ds = sintheta * math.cos(center_angle) - costheta * math.sin(center_angle)    
    dc = costheta * math.cos(center_angle) + sintheta * math.sin(center_angle)  
    dtheta = np.arctan2(ds,dc)

    orientation_component =  np.exp(-0.5 * (dtheta/angle_bandwidth)**2)

    ## frequency componenet ##
    # go to polar space
    raw = np.sqrt(x**2+y**2)
    # set origin to 1 as in the log space zero is not defined
    raw[mid,mid] = 1
    # go to log space
    raw = np.log2(raw)

    center_scale = math.log2(N) - f_0
    draw = raw-center_scale
    frequency_component = np.exp(-0.5 * (draw/ scale_bandwidth)**2)

    # reset origin to zero (not needed as it is already 0?)
    frequency_component[mid,mid] = 0

    return frequency_component * orientation_component

0
投票

是的 - 这对我也很有帮助。由于我刚刚学习 log-gabor 过滤器,我花了一些时间才意识到一些可能是基本的事情,但我想我应该在这里向其他人提及它。 Log-gabor 是频率滤波器。因此,您只需将其在频率空间中相乘,而不是制作一个小内核并在图像上进行卷积。

这意味着 N 必须与图像大小匹配,并且您必须首先使用 fft 进行变换,然后乘以滤波器,然后使用逆 fft 变换回来。这段代码为我提供了进行完整的端到端测试(过滤图像)所需的内容。然后我可以渲染生成的图像并验证我是否得到了我所期望的结果。作为测试,通过合并到单个 RGB 图像中可以轻松可视化 3 个方向。

def fft(im):
    return np.fft.fftshift(np.fft.fft2(im))

def ifft(f):
    return np.real(np.fft.ifft2(np.fft.ifftshift(f)))

N = src.shape[0]
f_0 = 4
number_orientations = 4
lg = [getLogGaborFilter(N,f_0,x,number_orientations) for x in range(0,number_orientations)]
f = [ifft(fft(src) * lg[x]) for x in range(0,number_orientations)]

我稍微更改了函数签名,以获取图像大小和方向数。另外,我花了一点时间才意识到 f_0 到底是什么。从 math.log2(N) - f_0 中,这相当于将 N 除以 2 的 f_0 次方,选择越来越大的尺度特征。我的理解是,这意味着您可能只想使用小的正整数。

© www.soinside.com 2019 - 2024. All rights reserved.