图像金字塔技术在很多层面上都有着广泛的应用，很多开源的工具也都有对他们的建立写了专门的函数，比如IPP，比如OpenCV等等，这方面的理论文章特别多，，但是大部多分开源的代码的实现都不是严格意义上的金字塔，而是做了一定的变通，这种变通常常为了快捷的实现类似的效果，本文将严格按照定义实现高斯金字塔的定义来优化该算法的实现。
图像金字塔技术在很多层面上都有着广泛的应用，很多开源的工具也都有对他们的建立写了专门的函数，比如IPP，比如OpenCV等等，这方面的理论文章特别多，我不需要赘述，但是我发现大部多分开源的代码的实现都不是严格意义上的金字塔，而是做了一定的变通，这种变通常常为了快捷的实现类似的效果，虽然这种变通不太会影响金字塔的效果，但是我这里希望从严格意义上对该算法进行优化，比如简要贴一下下面的某个高斯金字塔的代码：
上面的过程是对原图线进行双线性的下采样插值，然后再对采样后的图进行一定的高斯模糊。我们说这样做未尝不可，而真正的高斯金字塔的建立过程是：对上一级的数据进行高斯模糊（5x5)得到结果T，然后删除T的奇数行和奇数列数据作为下一级的结果（以0为下标起点的行列）。在此过程中使用到的高斯模糊的权重矩阵的形式如下所示：
应该说，上面的过程用伪代码表示应该是:
从编程优化的角度考虑，我们没有必要完整的对上一级进行高斯模糊，而只需要进行偶数行和偶数列的计算，然后赋值到下一层数据中，而进一步考虑上述矩阵的特殊性，可以通过减少重复计算以及合并相同系数项等手段来优化。对于边缘（2行2列）的像素，把他单独提取出来，减少中间数据的边缘判断可进一步提高速度。
再提出第一版C语言代码之前，我们来要看下各层金字塔的大小问题，如果上一层的大小位偶数，则下一层直接就除以2，但是如果是奇数，则除2时需要向上取整，比如宽某层的宽度尺寸为101，则之后的各层依次为101->51->26->13->7->4->2->1。
再看看前面权重矩阵的式子最右边那个乘法，那表示这个权重矩阵是行列可分离的，我们可以先计算行的加权，然后再利用这个加权值计算列的加权，也可以先计算列然后再计算行，这样原本每个像素处的25个乘法和多24次加法就可以减少为10次乘法和9次加法。对于沿着宽度方向的更新计算，我们还可以充分利用列方向的重叠累计信息，减少计算量，一个可行的简单的代码如下所示：
注意到，在单通道的示例里，我们用了5个中间变量，分别记录了某个位置5列像素的累加和，在移动到下一个目标像素时，由于我们是隔行隔列采样的，因此移动到下一个像素时只有3个位置时重叠的，也就是说只有3个分量可以重复利用，另外两个必须重新计算。
上面的代码是先计算列方向，然后在计算行方向的，基本上不需要额外的内存，我们再来试下先计算行方向的累积值，再处理列方向的方案：
其中IM_DownSampleLine8U_C函数如下所示：
注意到上述代码在结构上和第一版本其实差不多，不过多了5行临时内存，在更新行权重的时候也是只需要更新2行的，而无需整体更新。
我们对这两种方案进行了速度测试，由于本身这个的执行速度就比较块，因此我们对算法进行了100次计算，对于第一级为1920*1080大小的灰度图，下一级的高斯金字塔大小为960*540像素，算法C1测试的结果为267ms，算法C2的执行速度约为256ms，这说明他们本质上不存在大的速度差异（这里的时间都不包括处理边缘像素的，后面还要讨论，并且高斯金字塔也是采用unsigned char类型数据来保存的）。
在某些场合，我们还需要加快这个过程的速度，因此我考虑使用SSE优化他，考虑以上两种实现方式，哪一种更有利于SSE的处理呢，由于第一种方式前后的依赖比较强，用SSE做不是不可以，但估计效率不会有提升，需要太多次数据重组了，而第二种方式的由中间数据计算最后的结果很明显可以使用SSE处理，即下面的这三行代码：
这里的Sum是16位的，LinePD是8位的。替换的代码如下所示：
这么个简单的改动，速度大概到了180ms。
有点麻烦的是IM_DownSampleLine8U_C这个函数的优化，其核心代码如下：
最简单的SSE处理方式是加载5次不同位置的Src值，然后将数据转换为16位，再进行加法和乘法计算。
一次性处理8个像素，需要多次加载内存，原以为速度会有问题，结果一测，速度居然飙升到40ms，单次只需要0.4ms了。真的很高兴。
但是和普通的C比较一下，似乎结果不对啊，仔细分析，原来是因为这个对Src取样计算时每次是隔一个点取一个样的，而上述代码侧连续采样的，那怎么办呢？
也很简单，我们使用_mm_loadu_si128一次性加载16个字节，然后每隔一个像素就置0，这样就相当于把剩下的8个像素的值直接变为了16位数据了，一举两得。如下所示：
我们上面用的and运算来将有关位置0，当然还可以使用 shuffle指令来得到同样的结果，速度大概也就稍微慢一点，大概再45ms。
实际上，在这里的由于权重有一些特殊性，比如有2个1，2个4，4还可以用移位实现，如果是一些其他不太有规律的权重，比如 3 7 9 7 3这种，我们实际上还有一种优化方式来处理，因为在SSE里还有一个_mm_maddubs_epi16这个指令，他可以一次性实现16个字节数*16个signed char，然后再两两相加保存到8个short类型中去，比如上面的代码也可以用下面的方式实现：
在本例中，上述代码执行100次要50几毫秒，比前面的慢了，但是这里的组合确实是蛮有味道的，各种数据的灵活应用也是值得参考的。我反而更欣赏这段代码。
以上谈及的均是单通道的算法，如果是BGR 3个通道或者BGRA 4个通道的图像数据，情况就会复杂一些，但是同样的道理，可以使用shuffle来调整位置，然后使用类似的方式处理。
我们来谈谈浮点版本的高斯金字塔，这个再很多情况下也有需求，毕竟有很多算法是再浮点上进行处理的，那浮点版本的普通C的代码其实和C语言是差不多的，只需要将有关数据类型改为浮点就可以了，那对于其核心的DownSampleLine函数，也是我们优化的关键和难点，由于SSE一次性只能加载4个浮点数，如果还是和刚才处理字节数据那样，隔一个数取一个数，那么利用SSE一次性只能处理2个像素，而我们通过下面的美好的优化方式，一次性就能处理4个像素了，而且代码也很优美，我很是喜欢。
具体的每句代码的意思根据我上面的注释应该很容易能弄懂的，我就不加解释了。
最后，我们来关注下边缘的处理，边缘部分由于取样时会超出图像边界，因此，需要做判断，一种合理的方式是采用镜像数据，此时可以保证权重一定是256，我做了一个简单的函数：
为了程序完整，我们最后再加上周边像素的处理，然而我们发现一个严重的问题，再没有处理四周的函数中，我们运行100次SSE的耗时大约是45ms,一旦加入边缘像素的处理，这个耗时我们发现75ms，而普通C语言版本里由原来的260ms变为290ms，我们可能感受不到大的区别，但SSE的优化后，边缘部分居然占用了40%的耗时，因此，此时边缘特殊像素的处理就成了核心的事情了。
一种可行的优化方式就是类似于我前面做的Sobel边缘检测时方式，先对数据进行扩展，然后对扩展后的数据进行处理，此时边缘部分的处理已经被包括到SSE里去了，我尚未实践此方案的可行性和速度效果，相信应该不成问题。
附本文相关工程代码供参考：https://files.cnblogs.com/files/Imageshop/GaussPyramid.rar
