[人工智能]OpenCV学习笔记18_图像处理中的金字塔是个啥？

1、图像缩放

图像的缩放是线性变换的过程，在仿射变换中讲过，实际上是坐标变换之后再反向映射插值的过程，可以使用仿射矩阵，也可以使用OpenCV自带的resize函数，要缩小图像，一般推荐使用CV_INETR_AREA(区域插值)来插值；若要放大图像，推荐使用CV_INTER_LINEAR（线性插值）。

2、尺度空间

①高斯尺度空间(构建高斯尺度金字塔用)：

高斯尺度空间就是不同高斯标准差的滤波结果所构成的空间；

图中：Gaussian Space的${\sigma }_0=\sigma$,${\sigma }_1=k\sigma$,${\sigma }_2=k^2\sigma$…逐渐增大，构成一段连续的图像空间，${\sigma }_0,{\sigma }_1,{\sigma }_2,{\sigma }_3...$称为高斯尺度空间

②高斯差分尺度空间(DoG-Difference of Gaussian)：

高斯函数的差分，是灰度图像增强和角点检测的一种方法：
高斯尺度函数为：

$$G\left(x,y,\sigma \right)=\frac{1}{2\pi {\sigma }^2}{e}^{?\frac{{\left(x?u\right)}^2+{\left(y?v\right)}^2}{2{\sigma }^2}}$$
其中$\sigma$为高斯函数的标准差;
DoG算子表示为：
$$D\left(x,y,{\sigma }_2?{\sigma }_1\right)=[G\left(x,y,{\sigma }_2\right)?G\left(x,y,{\sigma }_1\right)]?f(x,y)$$
大意为：
标准差为${\sigma }_2$的高斯滤波器与图像的卷积结果减去标准差为${\sigma }_1$的图像的卷积结果就是Difference of Gaussian；
(根据卷积的结合率，可以得到卷积结果等于标准差为${\sigma }_2$的高斯滤波器减去标准差为${\sigma }_1$的滤波器再与原始图像卷积)
尺度空间理论目的是模拟图像数据的多尺度特征
尺度空间就是在分辨率相同的情况下，对于不同的高斯标准差构建不同尺度的高斯模糊图像

(1)DoG算子近似LoG算子：

.
LoG算子：
Laplace算子是通过对图像进行微分操作实现边缘检测的，所以对离散点和噪声比较敏感。于是，首先对图像进行高斯卷积滤波进行降噪处理，再采用Laplace算子进行边缘检测，就可以提高算子对噪声和离散点的鲁棒性，如此，拉普拉斯高斯算子Log（Laplace of Gaussian）就诞生了。

但是由于LoG算子在工程实现上计算量较大，所以一般采用逼近函数来近似。。。

假设${▽}^2$为拉普拉斯运算：
$$LoG={▽}^2{G}_{\sigma }(x,y,\sigma )?f(x,y)$$
其中$f(x,y)$表示原始图像

DoG算子：
$$D\left(x,y,{\sigma }_2?{\sigma }_1\right)=[G\left(x,y,{\sigma }_2\right)?G\left(x,y,{\sigma }_1\right)]?f(x,y)$$
则近似项为：
$$[G\left(x,y,{\sigma }_2\right)?G\left(x,y,{\sigma }_1\right)]\approx {▽}^2{G}_{\sigma }(x,y,\sigma )$$
也就是尺度差高斯核结果与二维高斯的拉普拉斯(二阶偏导)结果近似。

(2)推导：

三维空间的热传导方程的拉普拉斯形式为
$$\frac{?u}{?t}=k{▽}^{2}u$$
其中${▽}^{2}u$定义为u的拉普拉斯算子，$t$为时间，$u$为温度函数$u\left(x,y,z\right)$;
其中：
$${▽}^{2}u=\frac{{?}^{2}u}{?x^2}+\frac{{?}^{2}u}{?y^2}+\frac{{?}^{2}u}{?z^2}$$
假设高斯函数为$$G_{\sigma }\left(x,y,z,\sigma \right)$$
z是假设的参数,并不存在，只是为了对应方便;

则：
$$\frac{?G}{?\sigma }\approx k{▽}^{2}G$$
那么
$$\frac{?G}{?\sigma }=\frac{G\left(x,y,z,m\sigma \right)?G\left(x,y,z,\sigma \right)}{m\sigma ?\sigma }\approx k{▽}^{2}G$$
可知：
$$G\left(x,y,z,k\sigma \right)?G\left(x,y,z,\sigma \right)\approx (m?1)k?\sigma {▽}^{2}G$$
图像为：

由于我们只关心尺度空间的极值点，所以：
$$DoG\approx LoG$$

③分辨率空间：

分辨率空间是模拟人眼看物体的远近问题，人眼离物体越远物体越小，细节越模糊，人眼离物体越近物体越大，细节越明显。
例子：
可以想象一下你开车的时候前面有一辆黑色的车：
距离很远的时候你看见的是一团黑色，甚至分不清是什么；
稍微近一些你可以看出来这是一辆黑色的车，但是分不清具体是轿车还是suv；
再近一些你可以看出来这是一辆黑色的轿车，但是看不清楚是哪一款；
再近一些你可以看出来这是一辆奥迪RS7西装暴徒，轮廓优美；
再近一些。。。你就追尾了。。。。(车主会跟你详细介绍细节----分辨率越大，细节越明显)

例如：
530x390分辨率的图像为：

62X48分辨率的图像为：

在图像处理中分辨率可以通过上下采样获得插值后的分辨率，也可以通过坐标变换然后插值的方法获得等等；

3、图像金字塔(高斯金字塔)

百度百科的解释为：
图像金字塔是图像多尺度表达的一种，是一种以多分辨率来解释图像的有效但概念简单的结构。一幅图像的图像金字塔是一系列以金字塔形状（自下而上）逐步降低，且来源于同一张原始图的图像分辨率集合。其通过梯次向下采样获得，直到达到某个终止条件才停止采样。我们将一层一层的图像比喻成金字塔，层级越高，则图像越小，分辨率越低。
一幅图像的金字塔是一系列以金字塔形状排列的分辨率逐步降低，且来源同一张原始图的图像集合。
金字塔的底部是原始图像，而顶部是低分辨率的近似。
如下图所示：

具体操作为：
1、原始图像高斯滤波得到高斯模糊图像
2、高斯模糊图像下采样得到下一层图像
3、重复操作

结果为：

4、拉普拉斯金字塔

拉普拉斯金字塔称为图像金字塔的逆操作。
原始图像减去经过图像金字塔下采样之后再上采样的图像就是一层拉普拉斯操作；

具体操作为：
1、原始图像高斯滤波得到高斯模糊图像
2、高斯模糊图像下采样降低分辨率得到下一层图像
3、下一层图像上采样扩大分辨率得到上采样图像
4、上采样图像高斯滤波之后得到近似原始图像
5、原始图像减去上采样经过高斯滤波的图像得到一层拉普拉斯图像

结果：

5、下采样

图像的下采样为金字塔向上取图操作，金字塔越往上图像的分辨率越低，图像细节越少，下采样会丢失像素信息。。。
在OpenCV中的操作为：

CV_EXPORTS_W void pyrDown( 
InputArray src, OutputArray dst,
const Size& dstsize = Size(), 
int borderType = BORDER_DEFAULT );

参数简单介绍一下：
输入输出不用说；
dstsize为输出大小；
borderType 为边缘处理方法；

原始图像为：

原始图像下采样结果：

6、上采样

与金字塔方向相反，图像的上采样为金字塔图像向下取图操作，金字塔从下往上分辨率逐渐降低，上采样需要重组分辨率，也就是插值，是图像分辨率扩大的过程。。。

在OpenCV中的操作为：

CV_EXPORTS_W void pyrUp( 
InputArray src, OutputArray dst,
const Size& dstsize = Size(), 
int borderType = BORDER_DEFAULT );

参数跟下采样一样：
输入输出不用说；
dstsize为输出大小；
borderType 为边缘处理方法；

原始图像下采样后再上采样结果与原始图像对比：

7、高斯尺度金字塔(Sift特征常用)

高斯尺度金字塔分为三步:

1、构建高斯尺度空间；

(假设尺度空间为五组)
每组高斯空间的尺度值${\sigma }_0=\sigma$,${\sigma }_1=k\sigma$,${\sigma }_2=k^2\sigma$，${\sigma }_3=k^3\sigma$，${\sigma }_4=k^4\sigma$逐渐增大；

论文中取k值为$$k=2^{\frac{1}{s}}$$其中$$s=2$$(可能五组的时候效果比较好，总的组数=s+3)

所以k的初始值为$2^{\frac{1}{2}}$

2、构建分辨率空间；

步骤：
1、将第一步的高斯尺度空间依次在原始分辨率的图像上滤波作为第一层octave_0；

2、将分辨率降低为上一层图像的$\frac{1}{4}$作为第二层的分辨率，此时像素点行与列各降一半，使得二维高斯函数的标准差“作用域”扩大了4倍，所以需要二倍$\sigma$的初始值$2\sigma$才能使得尺度相同 。

行为$2\sigma$，列也为$2\sigma$，所以需要与上一层高斯组的倒数第三组图像尺度${\sigma }_2=k^2\sigma$作为初始尺度。

(例如：高斯核为5x5的核，图像缩小4倍，像素丢失，原来5x5的核现在的作用域实际上作用在了原始的10x10的图像上，所以需要扩大标准差近似达到10x10的范围)

3、往下每一层重复上一层操作。

3、最后构建高斯尺度金字塔；

论文中取k值为$$k=2^{\frac{1}{s}}$$其中$$s=2$$

所以k的初始值为$2^{\frac{1}{2}}$

结果：5层octave，5组高斯尺度，25张图像，$\sigma$=1.6

8、高斯差分金字塔(DoG)

在同一层中，每一张图像都是高斯尺度金字塔上一组尺度图像减去下一组尺度图像得到的图像；

DifferenceofGaussian正是表示高斯不同尺度下的差值，在SIFT特征算法中是取极值的重要步骤。

将差分值扩大70倍之后的差分金字塔：

参考文章：
SIFT算法：DOG尺度空间生产
SIFT特征提取分析

代码(不清楚的地方可以问我，要是有错误的地方欢迎指正)

//图像金字塔
#include <opencv.hpp>
#include <features2d/features2d.hpp>
using namespace std;
using namespace cv;

/*
NM_pyr 代表图像金字塔；

Gs_pyr 代表高斯金字塔；

DoG_pyr 代表高斯差分金字塔；

Lp_pyr 代表拉普拉斯金字塔；

*/

//普通图像金字塔
void Normal_Pyramid(Mat &src, vector<Mat>& NM_pyr, int nOctaves);

//高斯金字塔
void Gaussian_Pyramid(Mat &src, vector<Mat>& Gs_pyr, int nOctaves, double sigma = 1.6, int nOctavesLayers = 2);

//DoG金字塔
void DifferenceofGaussian(Mat &src,vector<Mat>& Gs_pyr, vector<Mat>& DoG_pyr, int nOctavesLayers=2);

//拉普拉斯金字塔
void Laplacian_Pyramid(vector<Mat>& NM_pyr, vector<Mat>& Lp_pyr);


void main()
{
	Mat src = imread("test.jpg");
	Mat dst_up, dst_down;


	vector<Mat> Gs_pyr,NM_pyr, Lp_pyr, DoG_pyr;
	//图像金字塔(普通高斯金字塔)
	Normal_Pyramid(src, NM_pyr, 5);
	//拉普拉斯金字塔
	Laplacian_Pyramid(NM_pyr, Lp_pyr);
	//高斯尺度金字塔
	Gaussian_Pyramid(src, Gs_pyr, 5);
	//DoG高斯差分金字塔
	DifferenceofGaussian(src, Gs_pyr, DoG_pyr);

	for (size_t i = 0; i < 5; i++)
	{
		cv::String str_1 = "NM_pyr";
		str_1 += to_string(i % 5);
		imshow(str_1, NM_pyr[i]);


	}
	for (size_t i = 0; i < 4; i++)
	{
		cv::String str_1 = "Lp_pyr";
		str_1 += to_string(i % 4);
		imshow(str_1, Lp_pyr[i]);


	}

	for (size_t i = 0; i < 25; i++)
	{
		cv::String str_1= "Gs_pyr";
		cv::String str_2 = "_";
		cv::String str;
		str_1 += to_string(i/5);
		str_2 += to_string(i%5);
		str = str_1 + str_2;
		imshow(str,Gs_pyr[i]);


	}
	for (size_t i = 0; i < 21; i++)
	{
		cv::String str_1 = "DoG_pyr";
		cv::String str_2 = "_";
		cv::String str;
		
		if (i < 5)
		{
			str_1 += to_string(i / 5);
			str_2 += to_string(i);
			str = str_1 + str_2;
			imshow(str, DoG_pyr[i]);

		}
		else
		{		
			str_1 += to_string(((i-1) / 4));
			str_2 += to_string((i-1)%4);
			str = str_1 + str_2;
			imshow(str, DoG_pyr[i]);
		}
		
		


	}
	pyrDown(src, dst_down);
	pyrUp(dst_down, dst_up);
	
	imshow("src", src);
	imshow("dst_up", dst_up);
	imshow("dst_down", dst_down);

	waitKey(0);
}

//普通图像金字塔
void Normal_Pyramid(Mat & src, vector<Mat>& NM_pyr, int nOctaves)
{
	//默认5x5的滤波器
	Mat kernel_x = (Mat_<float>(1, 5) << 0.125/2, 0.5/2, 0.75/2, 0.5/2, 0.125/2);//权值归一化处理
	Mat kernel_y = (Mat_<float>(5, 1) << 0.125/2, 0.5/2, 0.75/2, 0.5/2, 0.125/2);

	NM_pyr.resize(nOctaves);

	Mat dst_x,dst_xy,dst_uc,dst;
	for (size_t i = 0; i < nOctaves; i++)
	{
		if (i == 0)
		{
			NM_pyr[i] = src.clone();
		}
		else
		{	
		src.convertTo(dst_uc, CV_32F,1./255);
		//1、高斯快速滤波
		filter2D(dst_uc, dst_x,-1, kernel_x);
		filter2D(dst_x, dst_xy, -1, kernel_y);
		//2、分辨率降低
		dst_xy.convertTo(dst_xy, CV_8U, 255.0);

		resize(dst_xy, dst, Size(src.cols / pow(2, i), src.rows / pow(2, i)),0,0,INTER_NEAREST);
		
		NM_pyr[i] = dst;
		}
	}


}




//构建高斯尺度金字塔
void Gaussian_Pyramid(Mat &src, vector<Mat>& Gs_pyr, int nOctaves, double sigma, int nOctavesLayers )
{
	//1、构建高斯尺度空间 Layers(图层集合)
	//sigma的初始值，原始图像，尺度空间的层数
	vector<double> sig(nOctavesLayers+3);//构建sigma数组，每个图像一个sigma
	Gs_pyr.resize(nOctaves*(nOctavesLayers+3));//重新分配高斯金字塔的空间一共nOctaves层，每个分辨层的尺度为5层

	sig[0] = sigma;//初始值为输入sigma

	//nOctavesLayers = 2

	double k = pow(2., 1. / nOctavesLayers);// k = 2^(1/2)
	for (int i = 1 ; i < nOctavesLayers + 3 ;i++)//尺度空间层数为5
	{
	//构建sigma空间，sigma数组,每个图像一个sigma
	//sqrt（k^2*sigma^2  - sigma^2）
	//\sigma_{total}^2 = \sigma_{i}^2 - \sigma_{i-1}^2
		double sig_pre = pow(k, i)*sigma;
		double sig_total = k* sig_pre;
		sig[i] = sqrt(sig_total*sig_total - sig_pre*sig_pre);
	// 		sig[i] = k*sigma;		
	// 		k *= k;
	}
	int size_x = src.rows;
	int size_y = src.cols;
	//2、构建高斯分辨率空间
	for (size_t n = 0; n < nOctaves; n++)//图像层数为nOctaves   每层图像有nOctavesLayers + 3 个尺度空间
	{

		for (int i = 0 ;  i < nOctavesLayers + 3 ; i++)
		{
			Mat dst;
			resize(src, dst, Size(size_y / (pow(2, n)), size_x / (pow(2, n))),0,0, INTER_AREA);
			//构建尺度空间
			//sigma的值
			double sigma_x = pow(nOctavesLayers,n)*sig[i];
			double sigma_y = pow(nOctavesLayers, n)*sig[i];
			GaussianBlur(dst, dst, Size(5, 5), sigma_x, sigma_y);
			Gs_pyr[(n*(nOctavesLayers + 3)) + i] = dst;

		}


	}






}

//DoG金字塔
void DifferenceofGaussian(Mat &src,vector<Mat>& Gs_pyr, vector<Mat>& DoG_pyr, int nOctavesLayers)
{
	int n = Gs_pyr.size()/( nOctavesLayers + 3);
	int nO = nOctavesLayers + 2;
	DoG_pyr.resize(n*nO + 1);
	Mat dst;
	
	for (size_t i = 0; i < n; i++)
	{

		for (size_t j = 0; j < nO; j++)
		{
			if (i == 0 && j == 0)
			{
				DoG_pyr[i] = src - Gs_pyr[i];
				DoG_pyr[i + 1] = (Gs_pyr[i] - Gs_pyr[i + 1])*70;//为了便于观察，将差分值扩大70倍
			}
			else
			{
				DoG_pyr[i*nO +j+1] = (Gs_pyr[i*(nO+1) + j] - Gs_pyr[i*(nO + 1) + j + 1]) * 70;
				/*dst = DoG_pyr[i*nO + j + 1];*/
			}
			

		}

	}

}


// 拉普拉斯金字塔
// L = Gi - PyrUp(PyrDown(Gi))
void Laplacian_Pyramid(vector<Mat>& NM_pyr, vector<Mat>& Lp_pyr)
{
	//1、默认5x5的滤波器
	Mat kernel_x = (Mat_<float>(1, 5) << 0.125 / 2, 0.5 / 2, 0.75 / 2, 0.5 / 2, 0.125 / 2);
	Mat kernel_y = (Mat_<float>(5, 1) << 0.125 / 2, 0.5 / 2, 0.75 / 2, 0.5 / 2, 0.125 / 2);

	//2、构建拉普拉斯金字塔
	int nOctave = NM_pyr.size();//高斯金字塔的层数
	Lp_pyr.resize(NM_pyr.size()-1);//拉普拉斯金字塔的总层数

	Mat down, up, down_up;
	Mat dst_x, dst_xy,dst;
	for (size_t i = 0; i < nOctave-1; i++)
	{
		down = NM_pyr[i + 1];
		up = NM_pyr[i];

		resize(down, down_up, up.size());
		filter2D(down_up, dst_x, -1, kernel_x);
		filter2D(dst_x, dst_xy, -1, kernel_y);

		cv::subtract(up, dst_xy, dst);
		Lp_pyr[i] = dst;

	}





}