双三次插值及优化

1、数学模型对于一个目的像素，其坐标通过反向变换得到的在原图中的浮点坐标为(i+u,j+v),其中i、j均为非负整数，u、v为0,1)区间的浮点数，双三次插值考虑一个浮点坐标(i+u,j+v)周围的16个邻点，目的像素值f(i+u,j+v)可由如下插值公式得到：f(i+u,j+v)=A*B*CA=S(u+1)S(u+0)S(u-1)S(u-2)厂f(i-1,j-1)f(i-1,j+0)B=If(i+0,j-1)f(i+0,j+0)If(i+1,j-1)f(i+1,j+0)Lf(i+2,j-1)f(i+2,j+0)f(i-1,j+1)f(i-1,j+2)nf(i+0,j+1)f(i+0,j+2)If

2、(i+1,j+1)f(i+1,j+2)If(i+2,j+1)f(i+2,j+2)厂S(v+1)nC=IS(v+0)IIS(v-1)ILS(v-2)厂1-2*Abs(x)人2+Abs(x)人3,0v=Abs(x)v1S(x)=4-8*Abs(x)+5*Abs(x)A2-Abs(x)A3,1=Abs(x)=2S(x)是对Sin(x*Pi)/x的逼近(Pi是圆周率n)，为插值核。计算流程获取16个点的坐标P1、P2P16由插值核计算公式S(x)分别计算出x、y方向的插值核向量Su、Sv进行矩阵运算，得到插值结果iTemp1=Su0*P1+Su1*P5+Su2*P9+Su3*P13iTemp2=Su0

3、*P2+Su1*P6+Su2*P10+Su3*P14iTemp3=Su0*P3+Su1*P7+Su2*P11+Su3*P15iTemp4=Su0*P4+Su1*P8+Su2*P12+Su3*P16iResult=Sv1*iTemp1+Sv2*iTemp2+Sv3*iTemp3+Sv4*iTemp4在得到插值结果图后，我们发现图像中有“毛刺”，因此对插值结果做了个后处理，即：设该点在原图中的像素值为pSrc,若abs(iResult-pSrc)大于某阈值，我们认为插值后的点可能污染原图，因此用原像素值pSrc代替。3.算法优化由于双三次插值计算一个点的坐标需要其周围16个点，更有多达20次的乘法

4、及15次的加法，计算量可以说是非常大，势必要进行优化。我们选择了Intel的SSE2优化技术，它只支持在P4及以上的机器。测试当前CPU是否支持SSE2，可由CPUID指令得到，代码为：BOOLg_bSSE2=FALSE;_asmmoveax,1;cpuid;testedx,0 x04000000;jzNotSupport;movg_bSSE2,1NotSupport:支持SSE2的CPU引入了8个128位的寄存器，这样一个寄存器中就可以存放4个点（RGB）,有利于并行计算。详细代码见Transform.cpp中函数Optimize_Bicubic。优化中遇到的问题：图像每个点由RGB通道组成

5、，由于1个SSE2寄存器有16个字节，这样读入4个像素点后，要浪费4个字节，同时要花费时间将数据对齐，即由BRGB|RGBR|GBRG|BRGB对齐成0RGB|0RGB|0RGB|0RGB;读16字节数据到寄存器时，由于图像地址不能保证是16字节对齐，因此需用更多时钟周期的MOVDQU指令（6个以上时钟周期）;如能使地址16字节对齐,则可用MOVDQA指令（1个时钟周期）;为了消除除法及浮点运算，对权值放大256倍，这样在计算插值核时，必须用2Bytes来表示1个系数，而图像数据都是lByte,这样在对齐做乘法时，要浪费一半的SSE2寄存器的空间，导致运算时间变长;而若降低插值核的精度，使其在

6、1Byte表示范围内时，运算的精度又大为下降;对各指令的周期以及若干行指令是否能够并行流水缺乏经验和认识。附：SSE2指令整理算术（Arithmetic）指令:ADDPD-PackedDouble-PrecisionFloating-PointAddSSE22个double对应相加ADDPDxmm0,xmm1/m128SSEADDPS-PackedSingle-PrecisionFloating-PointAdd4个float对应相加ADDPSxmm0,xmm1/m128ADDSD-ScalarDouble-PrecisionFloating-PointAdd1个double（低端）对应相加S

7、SE2ADDSDxmm0,xmm1/m64ADDSS-ScalarSingle-PrecisionFloating-PointAddSSE1个float（低端）对应相加ADDSSxmm0,xmm1/m32PADDB/PADDW/PADDD-PackedAddOpcodeInstructionDescriptionOFFC/rPADDBmm,mm/m64Addpackedbyteintegersfrommm/m64andmm.660FFC/rPADDBxmm1,xmm2/m128Addpackedbyteintegersfromxmm2/m128andxmm1.0FFD/rPADDWmm,mm/

8、m64Addpackedwordintegersfrommm/m64andmm.660FFD/rPADDWxmml,xmm2/ml28Addpackedwordintegersfromxmm2/m128andxmm1.0FFE/rPADDDmm,mm/m64Addpackeddoublewordintegersfrommm/m64andmm.660FFE/rPADDDxmm1,xmm2/m128Addpackeddoublewordintegersfromxmm2/m128andxmm1.PADDQ-PackedQuadwordAddOpcodeInstructionDescriptionOF

9、D4/rPADDQmml,mm2/m64Addquadwordintegermm2/m64tomm166OFD4hrPADDQxmm1,xmm2/m128Addpackedquadwordintegersxmm2/m128toxmm1PADDSB/PADDSW-PackedAddwithSaturationOpcodeInstructionDescription0FEC/rPADDSBmm,mm/m64Addpackedsignedbyteintegersfrommm/m64andmmandsaturatetheresults.660FEC/rPADDSBxmm1,xmm2/m128Addpa

10、ckedsignedbyteintegersfromxmm2/m128andxmm1saturatetheresults.0FED/rPADDSWmm,mm/m64Addpackedsignedwordintegersfrommm/m64andmmandsaturatetheresults.660FED/rPADDSWxmm1,xmm2/m128Addpackedsignedwordintegersfromxmm2/m128andxmm1andsaturatetheresults.PADDUSB/PADDUSW-PackedAddUnsignedwithSaturationOpcodeInst

11、ructionDescription0FDC/rPADDUSBmm,mm/m64Addpackedunsignedbyteintegersfrommm/m64andmmandsaturatetheresults.660FDC/rPADDUSBxmm1,xmm2/m128Addpackedunsignedbyteintegersfromxmm2/m128andxmm1saturatetheresults.0FDD/rPADDUSWmm,mm/m64Addpackedunsignedwordintegersfrommm/m64andmmandsaturatetheresults.660FDD/rP

12、ADDUSWxmm1,xmm2/m128Addpackedunsignedwordintegersfromxmm2/m128toxmm1andsaturatetheresults.PMADDWD-PackedMultiplyandAddOpcodeInstructionDescription0FF5/rPMADDWDmm,mm/m64Multiplythepackedwordsinmmbythepackedwordsinmm/m64.Addthe32-bitpairsofresultsandstoreinmmasdoublewordX7X6X5x+xjX2X1xoY7Y6Y5Y4Y3Y2Y1Y

13、OABX7-Y7|AB9：X1-Y1iABXO-YO)OOHOOHOOHOOHOOHOOHSUMfTBWP7.TBWP0)PSADBW-PackedSumofAbsoluteDifferencesOpcodeInstructionDescriptionOFF6/rPSADBWmm1,mm2/m64Absolutedifferenceofpackedunsignedbyteintegersfrommm2/m64andmm1;differencesarethensummedtoproduceanunsignedwordintegerresult.66OFF6/rPSADBWxmml,xmm2/ml

14、28Absolutedifferenceofpackedunsignedbyteintegersfromxmm2/m128andxmm1;the8lowdifferencesand8highdifferencesarethensummedseparatelytoproducetwowordintegerresults.PSUBB/PSUBW/PSUBD-PackedSubtractOpcodeInstructionDescription0FF8/rPSUBBmm,mm/m64Subtractpackedbyteintegersinmm/m64frompackedbyteintegersinmm

15、.660FF8/rPSUBBxmm1,xmm2/m128Subtractpackedbyteintegersinxmm2/m128frompackedbyteintegersinxmm1.0FF9/rPSUBWmm,Subtractpackedwordintegersinmm/m64frommm/m64packedwordintegersinmm.66OFF9hrPSUBWxmml,xmm2/ml28Subtractpackedwordintegersinxmm2/m128frompackedwordintegersinxmm1.OFFA/rPSUBDmm,mm/m64Subtractpack

16、eddoublewordintegersinmm/m64frompackeddoublewordintegersinmm.66OFFAhrPSUBDxmml,xmm2/m128Subtractpackeddoublewordintegersinxmm2/mem128frompackeddoublewordintegersinxmm1.PSUBQ-PackedSubtractQuadwordOpcodeInstructionDescription0FFB/rPSUBQmm1,mm2/m64Subtractquadwordintegerinmm1frommm2/m64.660FFB/rPSUBQx

17、mm1,xmm2/m128Subtractpackedquadwordintegersinxmm1fromxmm2/m128.PSUBSB/PSUBSW-PackedSubtractwithSaturationOpcodeInstructionDescription0FE8/rPSUBSBmm,mm/m64Subtractsignedpackedbytesinmm/m64fromsignedpackedbytesinmmandsaturateresults.660FE8/rPSUBSBxmm1,xmm2/m128Subtractpackedsignedbyteintegersinxmm2/m1

18、28frompackedsignedbyteintegersinxmm1andsaturateresults.0FE9/rPSUBSWmm,mm/m64Subtractsignedpackedwordsinmm/m64fromsignedpackedwordsinmmandsaturateresults.660FE9/rPSUBSWxmm1,xmm2/m128Subtractpackedsignedwordintegersinxmm2/m128frompackedsignedwordintegersinxmm1andsaturateresults.PSUBUSB/PSUBUSW-PackedS

19、ubtractUnsignedwithSaturationOpcodeInstructionDescription0FD8/rPSUBUSBmm,mm/m64Subtractunsignedpackedbytesinmm/m64fromunsignedpackedbytesinmmandsaturateresult.660FPSUBUSBxmm1,SubtractpackedunsignedbyteintegersinD8/rxmm2/ml28xmm2/m128frompackedunsignedbyteintegersinxmmlandsaturateresult.OFD9/rPSUBUSW

20、mm,mm/m64Subtractunsignedpackedwordsinmm/m64fromunsignedpackedwordsinmmandsaturateresult.66OFD9/rPSUBUSWxmml,xmm2/m128Subtractpackedunsignedwordintegersinxmm2/m128frompackedunsignedwordintegersinxmm1andsaturateresult.SUBPD-PackedDouble-PrecisionFloating-PointSubtractOpcodeInstructionDescription660F5

21、C/rSUBPDxmm1,xmm2/m128Subtractpackeddouble-precisionfloating-pointvaluesinxmm2/m128fromxmm1.SUBPS-PackedSingle-PrecisionFloating-PointSubtractOpcodeInstructionDescription0F5C/rSUBPSxmm1xmm2/m128Subtractpackedsingle-precisionfloating-pointvaluesinxmm2/memfromxmm1.SUBSD-ScalarDouble-PrecisionFloating-

22、PointSubtractOpcodeInstructionDescriptionF20F5C/rSUBSDxmm1,xmm2/m64Subtractsthelowdouble-precisionfloating-pointnumbersinxmm2/mem64fromxmm1.SUBSS-ScalarSingle-FPSubtractOpcodeInstructionDescriptionF0F5CSUBSSxmm1,xmm2/m32Subtractthelowersingle-precisionfloating-pointnumbersinxmm2/m32fromxmm1.PMULHUW-

23、PackedMultiplyHighUnsignedOpcodeInstructionDescription0FE4/rPMULHUWmm1,Multiplythepackedunsignedwordintegersinmm1mm2/m64registerandmm2/m64,andstorethehigh16bitsoftheresultsinmm1.66OFE4/rPMULHUWxmml,xmm2/ml28Multiplythepackedunsignedwordintegersinxmm1andxmm2/m128,andstorethehigh16bitsoftheresultsinxm

24、m1.SRCX3X2X1XODESTY3Y2Y1YO73nY3Z2=X2tY2Z1=X1*Y1ZO=XO*YODEST23(31阴Zq3l-1GZl31-16Zq31-16PMULHW-PackedMultiplyHighSignedOpcodeInstructionDescription0FE5/rPMULHWmm,mm/m64Multiplythepackedsignedwordintegersinmm1registerandmm2/m64,andstorethehigh16bitsoftheresultsinmm1.660FE5/rPMULHWxmm1,xmm2/m128Multiply

25、thepackedsignedwordintegersinxmm1andxmm2/m128,andstorethehigh16bitsoftheresultsinxmm1.PMULLW-PackedMultiplyLowSignedOpcode0FD5/rInstructionDescriptionPMULLWmm,mm/m64Multiplythepackedsignedwordintegersinmm1registerandmm2/m64,andstorethelow16bitsoftheresultsinmm1.660FD5/rPMULLWxmm1,xmm2/m128Multiplythepackedsignedwordintegersinxmm1andxmm2/m128,andstorethelow16bitsoftheresultsinxmm1.PMULUDQ-MultiplyDoublewordUnsignedOpcodeInstructionDescription0FF4/rPMULUDQmm1,Multiplyunsigneddoublewordintegerinmm1bymm2/m64unsigneddoublewordintegerinmm2/m64,andstorethequadwordresultinmm1.66OFF4/