SSE2实现HAAR小波变换（dwt2与idwt2）

wiki链接：http://en.wikipedia.org/wiki/Haar_wavelet

 

可用SSE2实现HAAR小波变换，达到实时，关于HAAR小波的介绍可参考以上维基链接

 参考MATLAB中dwt2与idwt2的函数原型，基于OpenCV的框架进行了汇编优化实现

HAAR小波也可用于图像的压缩，将CH，CV，CD中的分量值小于某一阈值则归为0，从而这三个矩阵将成为稀疏矩阵（Sparse Matrix），反变换后的图像的质量将取决于选择阈值的大小。

 

实际实现时，可选择使用浮点数（单精度或双精度）进行矩阵计算，使用整数计算能得到更快的速度但不能进行完整的压缩与解压缩。

 

附代码：

 

[cpp]  view plain  copy

1. inline void dwt2_row(__out double* ca0,  
2.                      __out double* ch0,  
3.                      __out double* cv0,  
4.                      __out double* cd0,  
5.                      __in unsigned char* row0,  
6.                      __in unsigned char* row1,  
7.                      __in int col)  
8. {  
9.     __asm  
10.     {  
11.         mov         eax_ptr, ca0;  
12.         mov         ebx_ptr, ch0;  
13.         mov         ecx_ptr, cv0;  
14.         mov         edx_ptr, cd0;  
15.         mov         esi_ptr, row0;  
16.         mov         edi_ptr, row1;  
17.         pxor        xmm3, xmm3;  
18.         movapd      xmm7, g_halfd;  
19.         sub         col, 4;  
20.         jl          loop_2;  
21. loop_4:  
22.         movd        xmm1, [esi_ptr];  
23.         movd        xmm5, [edi_ptr];  
24.         punpcklbw   xmm1, xmm3;  
25.         punpcklbw   xmm5, xmm3;  
26.         punpcklwd   xmm1, xmm3;  
27.         punpcklwd   xmm5, xmm3;  
28.         cvtdq2pd    xmm0, xmm1;  
29.         cvtdq2pd    xmm4, xmm5;  
30.         shufpd      xmm1, xmm1, 1;  
31.         shufpd      xmm5, xmm5, 1;  
32.         cvtdq2pd    xmm1, xmm1;  
33.         cvtdq2pd    xmm5, xmm5;  
34.         addpd       xmm4, xmm0;  
35.         addpd       xmm5, xmm1;  
36.         mulpd       xmm4, xmm7;  
37.         mulpd       xmm5, xmm7;  
38.         subpd       xmm0, xmm4;  
39.         subpd       xmm1, xmm5;  
40.         movapd      xmm6, xmm4;  
41.         movapd      xmm2, xmm0;  
42.         shufpd      xmm4, xmm5, 0;  
43.         shufpd      xmm6, xmm5, 3;  
44.         shufpd      xmm0, xmm1, 0;  
45.         shufpd      xmm2, xmm1, 3;  
46.         addpd       xmm6, xmm4;  
47.         addpd       xmm2, xmm0;  
48.         mulpd       xmm6, xmm7;  
49.         mulpd       xmm2, xmm7;  
50.         subpd       xmm4, xmm6;  
51.         subpd       xmm0, xmm2;  
52.         movupd      [eax_ptr], xmm6;  
53.         movupd      [ebx_ptr], xmm4;  
54.         movupd      [ecx_ptr], xmm2;  
55.         movupd      [edx_ptr], xmm0;  
56.         add         esi_ptr, 4;  
57.         add         edi_ptr, 4;  
58.         add         eax_ptr, 0x10;  
59.         add         ebx_ptr, 0x10;  
60.         add         ecx_ptr, 0x10;  
61.         add         edx_ptr, 0x10;  
62.         sub         col, 4;  
63.         jge         loop_4;  
64. loop_2:  
65.         cmp         col, -2;  
66.         jl          loop_end;  
67.         pinsrw      xmm0, [esi_ptr], 0;  
68.         pinsrw      xmm4, [edi_ptr], 0;  
69.         punpcklbw   xmm0, xmm3;  
70.         punpcklbw   xmm4, xmm3;  
71.         punpcklwd   xmm0, xmm3;  
72.         punpcklwd   xmm4, xmm3;  
73.         cvtdq2pd    xmm0, xmm0;  
74.         cvtdq2pd    xmm4, xmm4;  
75.         addpd       xmm4, xmm0;  
76.         mulpd       xmm4, xmm7;  
77.         subpd       xmm0, xmm4;  
78.         movapd      xmm5, xmm4;  
79.         shufpd      xmm4, xmm0, 0;  
80.         shufpd      xmm5, xmm0, 3;  
81.         addpd       xmm5, xmm4;  
82.         mulpd       xmm5, xmm7;  
83.         subpd       xmm4, xmm5;  
84.         movsd       [eax_ptr], xmm5;  
85.         shufpd      xmm5, xmm5, 1;  
86.         movsd       [ebx_ptr], xmm4;  
87.         shufpd      xmm4, xmm4, 1;  
88.         movsd       [ecx_ptr], xmm5;  
89.         movsd       [edx_ptr], xmm4;  
90. loop_end:  
91.     }  
92. }  
93.   
94. inline void idwt2_row(__out unsigned char* row0,  
95.                       __out unsigned char* row1,  
96.                       __in double* ca0,  
97.                       __in double* ch0,  
98.                       __in double* cv0,  
99.                       __in double* cd0,  
100.                       __in int col)  
101. {  
102.     __asm  
103.     {  
104.         mov         eax_ptr, ca0;  
105.         mov         ebx_ptr, ch0;  
106.         mov         ecx_ptr, cv0;  
107.         mov         edx_ptr, cd0;  
108.         mov         esi_ptr, row0;  
109.         mov         edi_ptr, row1;  
110.         sub         col, 4;  
111.         jl          loop_2;  
112. loop_4:  
113.         movupd      xmm0, [eax_ptr];  
114.         movupd      xmm1, [ebx_ptr];  
115.         movupd      xmm4, [ecx_ptr];  
116.         movupd      xmm5, [edx_ptr];  
117.         addpd       xmm1, xmm0;  
118.         addpd       xmm5, xmm4;  
119.         addpd       xmm0, xmm0;  
120.         addpd       xmm4, xmm4;  
121.         subpd       xmm0, xmm1;  
122.         subpd       xmm4, xmm5;  
123.         movapd      xmm2, xmm1;  
124.         movapd      xmm6, xmm5;  
125.         shufpd      xmm1, xmm0, 0;  
126.         shufpd      xmm2, xmm0, 3;  
127.         shufpd      xmm5, xmm4, 0;  
128.         shufpd      xmm6, xmm4, 3;  
129.         addpd       xmm5, xmm1;  
130.         addpd       xmm6, xmm2;  
131.         addpd       xmm1, xmm1;  
132.         addpd       xmm2, xmm2;  
133.         subpd       xmm1, xmm5;  
134.         subpd       xmm2, xmm6;  
135.         cvttpd2dq   xmm5, xmm5;  
136.         cvttpd2dq   xmm6, xmm6;  
137.         cvttpd2dq   xmm1, xmm1;  
138.         cvttpd2dq   xmm2, xmm2;  
139.         shufpd      xmm5, xmm6, 0;  
140.         shufpd      xmm1, xmm2, 0;  
141.         packssdw    xmm5, xmm1;  
142.         packuswb    xmm5, xmm5;  
143.         pshufd      xmm1, xmm5, 1;  
144.         movd        [esi_ptr], xmm5;  
145.         movd        [edi_ptr], xmm1;  
146.         add         esi_ptr, 4;  
147.         add         edi_ptr, 4;  
148.         add         eax_ptr, 0x10;  
149.         add         ebx_ptr, 0x10;  
150.         add         ecx_ptr, 0x10;  
151.         add         edx_ptr, 0x10;  
152.         sub         col, 4;  
153.         jge         loop_4;  
154. loop_2:  
155.         cmp         col, -2;  
156.         jl          loop_end;  
157.         movsd       xmm0, [eax_ptr];  
158.         movsd       xmm1, [ebx_ptr];  
159.         movsd       xmm4, [ecx_ptr];  
160.         movsd       xmm5, [edx_ptr];  
161.         addpd       xmm1, xmm0;  
162.         addpd       xmm5, xmm4;  
163.         addpd       xmm0, xmm0;  
164.         addpd       xmm4, xmm4;  
165.         subpd       xmm0, xmm1;  
166.         subpd       xmm4, xmm5;  
167.         shufpd      xmm1, xmm0, 0;  
168.         shufpd      xmm5, xmm4, 0;  
169.         addpd       xmm5, xmm1;  
170.         addpd       xmm1, xmm1;  
171.         subpd       xmm1, xmm5;  
172.         cvttpd2dq   xmm5, xmm5;  
173.         cvttpd2dq   xmm1, xmm1;  
174.         packssdw    xmm5, xmm1;  
175.         packuswb    xmm5, xmm5;  
176.         movd        eax_ptr, xmm5;  
177.         mov         [esi_ptr], ax;  
178.         shr         eax_ptr, 16;  
179.         stosw;  
180. loop_end:  
181.     }  
182. }  
183.   
184. inline void dwt2(__out cv::Mat& CA,  
185.                  __out cv::Mat& CH,  
186.                  __out cv::Mat& CV,  
187.                  __out cv::Mat& CD,  
188.                  __in cv::Mat const& I)  
189. {  
190.     if(CA.type() != CV_64FC1 || CH.type() != CV_64FC1 || CV.type() != CV_64FC1 || CD.type() != CV_64FC1 || I.channels() != 1)  
191.         return;  
192.   
193.     double* ca = reinterpret_cast<double*>(CA.data);  
194.     double* ch = reinterpret_cast<double*>(CH.data);  
195.     double* cv = reinterpret_cast<double*>(CV.data);  
196.     double* cd = reinterpret_cast<double*>(CD.data);  
197.     unsigned char* row = reinterpret_cast<unsigned char*>(I.data);  
198.   
199.     for(int i=0; i<I.rows; i+=2)  
200.     {  
201.         dwt2_row(ca, ch, cv, cd, row, row+I.cols, I.cols);  
202.         ca += CA.cols;  
203.         ch += CH.cols;  
204.         cv += CV.cols;  
205.         cd += CD.cols;  
206.         row += I.cols*2;  
207.     }  
208. }  
209.   
210. inline void idwt2(__out cv::Mat& I,  
211.                   __in cv::Mat const& CA,  
212.                   __in cv::Mat const& CH,  
213.                   __in cv::Mat const& CV,  
214.                   __in cv::Mat const& CD)  
215. {  
216.     if(CA.type() != CV_64FC1 || CH.type() != CV_64FC1 || CV.type() != CV_64FC1 || CD.type() != CV_64FC1 || I.channels() != 1)  
217.         return;  
218.   
219.     double* ca = reinterpret_cast<double*>(CA.data);  
220.     double* ch = reinterpret_cast<double*>(CH.data);  
221.     double* cv = reinterpret_cast<double*>(CV.data);  
222.     double* cd = reinterpret_cast<double*>(CD.data);  
223.     unsigned char* row = reinterpret_cast<unsigned char*>(I.data);  
224.   
225.     for(int i=0; i<I.rows; i+=2)  
226.     {  
227.         idwt2_row(row, row+I.cols, ca, ch, cv, cd, I.cols);  
228.         ca += CA.cols;  
229.         ch += CH.cols;  
230.         cv += CV.cols;  
231.         cd += CD.cols;  
232.         row += I.cols*2;  
233.     }  
234. }  

版权归作者所有，转载请注明出处！