矢量化双线性插值算法


双线性插值算法由于速度较快、效果尚可,图像大小变换中经常使用,很多人都需要在软件中针对不同格式的数据实现。网上有不少双线性插值算法的实现方法,也有矢量化的一些方法,但是都不是很深入,对数据格式依赖较重,并行性也不高。

由于双线性插值算法是行列可分离的,在实现中使用行列分离可降低计算复杂度,因此各种实现也都是行列分离的。由于双线性插值的特性,在列方向上要实现矢量化还是比较容易的,也有各种算法介绍。但是在行方向上,由于双线性插值读取源数据近似随机的特性,基本上没有介绍矢量化算法的文章,有的也是多通道图像的多个通道的矢量化。如RGBA格式的图像,在行方向上同时计算RGBA四个通道的值,可以一定程度上矢量化。但是,由于图像处理中通常使用的格式是通道分离的格式、而不是通道组合的格式,这样的算法可以应用的场景是比较有限的。目前OpenCV 3.1的实现中,双线性插值的行方向也没有使用矢量化算法,只在列方向使用了矢量化算法。

如果使用SSSE3中新加入的查表指令pshufb,可以把随机位置读取源数据的操作矢量化。利用该指令,可以把一定参数范围内的行方向的双线性插值算法矢量化。基本原理是,对连续的n个输出列,计算其所有值的第一源数据位置srcPos[0]~srcPos[n]。根据双线性插值的原理,srcPos[0]~srcPos[n]必是单调递增的(srcPos[i]>=srcPos[i-1])。因此,计算n列输出的源数据的位置差最大为posDiff=srcPos[n]-srcPos[0]。如果posDiff不超过15,则可以读取以srcPos[0]开始的16个像素,然后以srcPos[0]~srcPos[n]为索引执行pshufb指令,即可读取到计算n列输出所需的n个第一源数据。根据双线性插值原理,第二源数据的索引正好是第一元数据的索引加1,因此同样原理可以读取到第二元数据。如此获得的源数据是按照输出列顺序排列的,因此完全可以用矢量指令并行运算。该算法的并行度限制为n列输出所需的源数据列差异不超过15。一般来说,因为需要使用16位整数运算,使用SSE2指令的双线性插值算法通常同时计算8个像素。如果n=8,则在精心设计的算法上,只要插值比例wOld/wNew<15/7~=2.142857,算法都可以正确运行。这意味着,只要插值行方向是放大或压缩不要超过15/7,算法都可以实现一次计算8像素的并行。如果n取4即一次只计算4个像素,则可实现最大的压缩比例为15/3。如果是RGBA的多通道格式,那么一次运算8个数据的话只想要2列并行,可支持最大3/1的压缩比例,不是15/1的原因是16字节寄存器此时只能装载4个像素的数据。如果用256位AVX指令,则可以支持压缩比例不超过31/15的16像素并行计算,或31/7的8像素并行计算。

上面的算法隐含了一个前提没有提,就是srcPos[0]~srcPos[n]的计算也必须矢量化。对第I列,其源数据位置应为I*delta,其中delta是插值比例wOld/wNew。于是,srcPos[i]=I*delta-srcPos[0],这里i是列在当前组的下标,而I是列在整个目的图像的下标。最直接的做法是根据这个公式使用pmulhw计算。不过,这样计算相对复杂,更高效的做法是迭代计算srcPos[0]。后一组的srcPos[0]是前一组的srcPos[n]+delta,因此循环中可迭代更新srcPos[0],则srcPos[0]~srcPos[n]的增量部分可事先算好,循环中并行加上srcPos[0]的迭代值即可。而且,srcPos[0]中整数部分可以总是分离出来加到源数据指针上,只保留小数部分用于迭代。这样,迭代中很容易计算出源数据下标,且同时获得的小数部分正好是插值权重。当然,由于计算的整数部分要占用4位(最大值15),小数部分只剩下12位了,插值权重的精度会差一些。

使用以上算法还可以有一个优化。由于行列分离计算,获得最终结果的计算可以先行后列,也可以先列后行。由于行插值复杂度较大,因此我们总是希望减少行插值的次数。因此,如果列方向是放大,则应采用先行后列的顺序,使用两行宽度为目的图像宽的中间缓冲区以适当算法缓冲数据可保证行插值次数为源图像列数。如果列方向是缩小,应采用先列后行的顺序,此时行插值次数为目的图像列数,需要一行源图像宽度的中间缓冲区。注意应避免用全幅图像的中间缓冲区将行列计算完全分离。因为这样的做法中间缓冲区无法存放在cache中,而1~2行的中间缓冲区可以存放在cache中,可以提升速度。

实现以上算法在Core I5-3317U 1.7GHz上测试,单线程大图像(4k*4k)插值,可实现0.97~1.2G/s的像素填充率,只做列方向的插值可实现2G/s以上的像素填充率。OpenCV 3.1的实现由于没有行插值的矢量化、也没有对行列顺序进行优化,在同样行列比例不同时性能差异很大。行压缩、列放大时性能最差,只能达到约0.12G/s的像素填充率。行放大、列压缩时,OpenCV 3.1可实现高达0.75G/s的像素填充率。其它情况下OpenCV 3.1性能有0.3-0.6G/s像素填充率。不过,上面算法的中间缓冲数据使用的8位,精度较差,OpenCV的实现使用了更高位数的中间数据,精度较高。由于我暂不需要这一点精度,这里的实现还是选取的8位中间缓冲区。速度测试中,二者都把中间缓冲区和目的图像的内存分配的时间计算在内,因为OpenCV的实现无法剥离这二者。

理论上,以上算法实现两次立方插值的矢量化也是可行的,但因为我暂不需要,所以这里没有实现。以上算法的源码可在这里下载(提取码h533)。汇编函数在ASM目录,是nasm格式,需要nasmx-1.4。C++代码在根目录下。由于各人环境差异较大,因此没有配置编译环境,自己加入你的工程即可。汇编原则上可支持Win32/Win64/Linux,但只在win32/win64测试过,Linux尚未测试。