OpenCV：拉普拉斯金字塔+Gabor方向滤波器原理及实现

一、拉普拉斯金字塔

1.1 原理

        拉普拉斯金字塔将图像分解为低频分量和多个不同尺度的高频带通分量。使用拉普拉斯金字塔可以提取到图像不同尺度的空间频率特征。

        拉普拉斯金字塔滤波器组首先将图像下采样，从而获取该尺度的低频分量。然后，将下采样后的图像插值进行上采样，将插值后的图像与原图像做差，得到的便是原图像尺度下的高频残差，也就是高频分量。对于低频分量不断地迭代这一过程，从而得到不同尺度的高频分量（带通）和低频分量。

二、Gabor滤波器

2.1 原理

        在图像处理中，Gabor函数是一个用于边缘提取的线性滤波器。Gabor滤波器的频率和方向表达同人类视觉系统类似。研究发现，Gabor滤波器十分适合纹理表达和分离。在空间域中，一个二维Gabor滤波器是一个由正弦平面波调制的高斯核函数。Gabor滤波器的表达式如下：
复数表达：

实数部分：

虚数部分：

其中：

下面介绍公式中各个参数的含义：

波长（λ）：它的值以像素为单位指定，通常大于等于2.但不能大于输入图像尺寸的五分之一。

方向（θ）：这个参数指定了Gabor函数并行条纹的方向，它的取值为0到360度

相位偏移（φ）：它的取值范围为-180度到180度。其中，0he180度分别对应中心对称的center-on函数和center-off函数，而-90度和90度对应反对称函数。

长宽比（γ）：空间纵横比，决定了Gabor函数形状（support，我翻译为形状）的椭圆率（ellipticity）。当γ= 1时，形状是圆的。当γ< 1时，形状随着平行条纹方向而拉长。通常该值为0.5

带宽（b）：Gabor滤波器的半响应空间频率带宽b和σ/ λ的比率有关，其中σ表示Gabor函数的高斯因子的标准差。

σ的值不能直接设置，它仅随着带宽b变化。带宽值必须是正实数，通常为1，此时，标准差和波长的关系为：σ= 0.56 λ。带宽越小，标准差越大，Gabor形状越大，可见平行兴奋和抑制区条纹数量越多。

Gabor滤波器实数部分如下所示：

Gabor滤波器虚数部分如下所示：

三、图像多尺度多方向分解

        使用拉普拉斯金字塔结合Gabor方向滤波器可以将图像分解到多尺度和多方向，代码实现如下：

#include <opencv2/opencv.hpp>
#include <iostream>

using namespace cv;
using namespace std;

void pyramid_up(Mat &image, vector<Mat> &pyramid_images, int level);
void laplaian_demo(vector<Mat> &pyramid_images, Mat &image, vector<Mat> &pyramid_gabor_images);
//int main(int artc, char** argv) {
int LaplaianTest(){
	Mat src = imread("D:\\...test.jpg");
	if (src.empty()) {
		printf("could not load image...\n");
		return -1;
	}
	namedWindow("input", WINDOW_AUTOSIZE);
	imshow("input", src);

	vector<Mat> p_images;
	vector<Mat> pyramid_gabor_images;
	pyramid_up(src, p_images, 3);
	laplaian_demo(p_images, src,pyramid_gabor_images);

	waitKey(0);
	return 0;
}

void pyramid_up(Mat &image, vector<Mat> &pyramid_images, int level) {
	Mat temp = image.clone();
	Mat dst;
	for (int i = 0; i < level; i++) {
		pyrDown(temp, dst);
		//imshow(format("pyramid_up_%d", i), dst);
		temp = dst.clone();
		pyramid_images.push_back(temp);
	}
}

void laplaian_demo(vector<Mat> &pyramid_images, Mat &image, vector<Mat> &pyramid_gabor_images) {
	double theta[4];

	// theta 法线方向
	theta[0] = 0;
	theta[1] = CV_PI / 4;
	theta[2] = CV_PI / 2;
	theta[3] = CV_PI - CV_PI / 4;

	int kernel_size = 3;
	double sigma = 1.0, lambd = CV_PI / 8, gamma = 0.5, psi = 0;
	Mat kernel[4];
	for (int ang = 0; ang < 4; ang++) {
		kernel[ang] = getGaborKernel(cv::Size(kernel_size, kernel_size), sigma, theta[ang], lambd, gamma, psi, CV_32F);
	}

	
	Mat temp;
	for (int t = pyramid_images.size() - 1; t > -1; t--) {
		Mat dst;
		if (t - 1 < 0) {
			pyrUp(pyramid_images[t], dst, image.size());
			subtract(image, dst, dst);
			dst = dst + Scalar(127, 127, 127);
			//imshow(format("laplaian_layer_%d", t), dst);
		}
		else {
			pyrUp(pyramid_images[t], dst, pyramid_images[t - 1].size());
			subtract(pyramid_images[t - 1], dst, dst);
			dst = dst + Scalar(127, 127, 127);
			//imshow(format("laplaian_layer_%d", t), dst);
		}
		
		for (int i = 0; i < 4; i++) {
			filter2D(dst, temp, CV_32F, kernel[i]);
			pyramid_gabor_images.push_back(temp);
			imshow(format("laplaian_layer_%d_d_%d", t, i), temp);
		}
		
	}
}

参考：Gabor滤波器学习