Generating Any Levels Le Gall 5/3 Cdf 5/3 Wavelet Transform Matrix For 1D Signal With Whole Point Symmetric Padding By Matlab

Le Gall 5/3 also named Cdf 5/3 Wavelet Transform is used by JPEG2k standard.

This is the matlab code that would generate any levels Le Gall 5/3 Cdf 5/3 Wavelet Transform Matrix for 1D Signal With Whole Point Symmetric Padding:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
function H = legall53td(X,nlevel)
 
% Author: Jun Li, more info@ http://goldensectiontransform.org/
% Generating any levels Le Gall 5/3 CDF 5/3 wavelet transform matrix
% for 1D signal with whole point symmetric padding by matlab
% X is of length 2^N and note that CDF 5/3 matrix is biorthogonal.
 
[xx,yy] = size(X);
 
H=X;
for i=1:nlevel
 
   for k=1:yy
      [ss,dd] = legall531d(H(1:xx,k)'); % column transform
      H(1:xx,k) = [ss,dd]';
   end
 
   xx = xx/2;
   yy = yy;
 
end
 
 
%% 1d LeGall 5/3 lifting scheme % symw_ext->...12321...
 
function [ss,dd] = legall531d(S)
 
N = length(S);
 
ga = -1/2;
gb = 1/4;
gc = sqrt(2);
 
s0 = S(1:2:N-1);  % S(1),S(3),S(5),S(7)...
d0 = S(2:2:N);    % S(2),S(4),S(6),S(8)...
 
d1 = d0 + ga*(s0 + [s0(2:N/2) s0(N/2)]);
s1 = s0 + gb*(d1 + [d1(1) d1(1:N/2-1)]);
 
ss = gc*s1;
dd = d1/gc;
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
function X = ilegall53td(H,nlevel)
 
% Author: Jun Li, more info@ http://goldensectiontransform.org/
% Generating any levels inverse Le Gall 5/3 CDF 5/3 wavelet transform matrix
% for 1D signal with whole point symmetric padding by matlab
% H is of length 2^N and note that CDF 5/3 matrix is biorthogonal.
 
[x,y] = size(H);
xx = x/(2^nlevel);
yy = y/(2^nlevel);
 
X=H;
for i=1:nlevel
 
   for j=1:yy*2
      ss = X(1:xx,j);
      dd = X(xx+1:xx*2,j);
      X(1:xx*2,j) = ilegall531d(ss',dd')'; % inverse column transform
   end
 
   xx = xx*2;
   yy = yy*2;
 
end
 
 
%% inverse 1d LeGall 5/3 lifting scheme % symw_ext->...12321...
 
function S = ilegall531d(ss,dd)
 
N = length(ss) + length(dd);
 
ga = -1/2;
gb = 1/4;
gc = sqrt(2);
 
d1 = gc*dd;
s1 = ss/gc;
 
s0 = s1 - gb*(d1 + [d1(1) d1(1:N/2-1)]);
d0 = d1 - ga*(s0 + [s0(2:N/2) s0(N/2)]);
 
S(2:2:N) = d0;
S(1:2:N-1) = s0;
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
function outputcdf53
 
% wh = 2;
% E = eye(wh);
% nlevel = 1;
 
wh = 8;
E = eye(wh);
nlevel = 3;
 
LG = legall53td(E,nlevel)
 
ILG = ilegall53td(E,nlevel)
 
%% output transform matrix to txt file
 
fid=fopen('C:\legall53mtx.txt','w');
 
fprintf(fid,'%%%% %g-level %g point CDF 5/3 matrix\r\n\r\n',nlevel,wh);
 
%% output CDF 5/3 wavelet transform matrix
 
fprintf(fid,'legall53mtx_%g_%g = [',wh,nlevel);
matrix=LG;
[m,n]=size(matrix);
for i=1:1:m
   for j=1:1:n
      if j==n
         fprintf(fid,'%g\r\n',matrix(i,j));
      else
         fprintf(fid,'%g ',matrix(i,j));
      end
   end
end
fprintf(fid,'];\r\n\r\n\r\n');
%% end
 
 
%% output inverse CDF 5/3 wavelet transform matrix
 
fprintf(fid,'%%%% %g-level %g point inverse CDF 5/3 matrix\r\n\r\n',nlevel,wh);
 
fprintf(fid,'ilegall53mtx_%g_%g = [',wh,nlevel);
matrix=ILG;
[m,n]=size(matrix);
for i=1:1:m
   for j=1:1:n
      if j==n
         fprintf(fid,'%g\r\n',matrix(i,j));
      else
         fprintf(fid,'%g ',matrix(i,j));
      end
   end
end
fprintf(fid,'];\r\n');
 
fclose(fid);
 
% is LG*ILG equal identity matrix:
 
isequal(round(legall53td(E,nlevel)*ilegall53td(E,nlevel)),E)

Results:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
%% 2 point
legall53mtx_2_1 = [sqrt(1/2) sqrt(1/2)
-sqrt(1/2) sqrt(1/2)
];
 
 
ilegall53mtx_2_1 = [sqrt(1/2) -sqrt(1/2)
sqrt(1/2) sqrt(1/2)
];
 
%% 4 point
legall53mtx_4_2 = [5/8 3/4 3/8 1/4
-7/8 -1/4 7/8 1/4
-sqrt(1/8) sqrt(1/2) -sqrt(1/8) 0
0 0 -sqrt(1/2) sqrt(1/2)
];
 
 
ilegall53mtx_4_2 = [1/2 -1/2 -sqrt(1/2) 0
1/2 0 sqrt(25/32) -sqrt(1/32)
1/2 1/2 -sqrt(1/8) -sqrt(1/8)
1/2 1/2 -sqrt(1/8) sqrt(9/8)
];
 
%% 8 point 1 level
legall53mtx_8_1 = [sqrt(9/8) sqrt(1/2) -sqrt(1/8) 0 0 0 0 0
-sqrt(1/32) sqrt(1/8) sqrt(9/8) sqrt(1/8) -sqrt(1/32) 0 0 0
0 0 -sqrt(1/32) sqrt(1/8) sqrt(9/8) sqrt(1/8) -sqrt(1/32) 0
0 0 0 0 -sqrt(1/32) sqrt(1/8) sqrt(25/32) sqrt(1/8)
-sqrt(1/8) sqrt(1/2) -sqrt(1/8) 0 0 0 0 0
0 0 -sqrt(1/8) sqrt(1/2) -sqrt(1/8) 0 0 0
0 0 0 0 -sqrt(1/8) sqrt(1/2) -sqrt(1/8) 0
0 0 0 0 0 0 -sqrt(1/2) sqrt(1/2)
];
 
 
ilegall53mtx_8_1 = [sqrt(1/2) 0 0 0 -sqrt(1/2) 0 0 0
sqrt(1/8) sqrt(1/8) 0 0 sqrt(25/32) -sqrt(1/32) 0 0
0 sqrt(1/2) 0 0 -sqrt(1/8) -sqrt(1/8) 0 0
0 sqrt(1/8) sqrt(1/8) 0 -sqrt(1/32) sqrt(9/8) -sqrt(1/32) 0
0 0 sqrt(1/2) 0 0 -sqrt(1/8) -sqrt(1/8) 0
0 0 sqrt(1/8) sqrt(1/8) 0 -sqrt(1/32) sqrt(9/8) -sqrt(1/32)
0 0 0 sqrt(1/2) 0 0 -sqrt(1/8) -sqrt(1/8)
0 0 0 sqrt(1/2) 0 0 -sqrt(1/8) sqrt(9/8)
];
 
 
%% 8 point 3 level
legall53mtx_8_3 = [sqrt(9/32) sqrt(1/2) sqrt(529/2048) sqrt(81/512) sqrt(25/512) sqrt(25/512) sqrt(49/2048) sqrt(1/128)
-sqrt(25/32) -sqrt(1/2) -sqrt(25/2048) sqrt(25/512) sqrt(441/512) sqrt(81/512) sqrt(9/2048) sqrt(1/128)
-1/2 0 15/16 1/8 -1/2 -1/8 1/16 0
0 0 1/8 -1/4 -7/8 0 3/4 1/4
-sqrt(1/8) sqrt(1/2) -sqrt(1/8) 0 0 0 0 0
0 0 -sqrt(1/8) sqrt(1/2) -sqrt(1/8) 0 0 0
0 0 0 0 -sqrt(1/8) sqrt(1/2) -sqrt(1/8) 0
0 0 0 0 0 0 -sqrt(1/2) sqrt(1/2)
];
 
 
ilegall53mtx_8_3 = [sqrt(1/8) -sqrt(1/8) -1/2 0 -sqrt(1/2) 0 0 0
sqrt(1/8) -sqrt(1/32) 1/16 -1/16 sqrt(25/32) -sqrt(1/32) 0 0
sqrt(1/8) 0 5/8 -1/8 -sqrt(1/8) -sqrt(1/8) 0 0
sqrt(1/8) sqrt(1/32) 3/16 -3/16 -sqrt(1/32) sqrt(9/8) -sqrt(1/32) 0
sqrt(1/8) sqrt(1/8) -1/4 -1/4 0 -sqrt(1/8) -sqrt(1/8) 0
sqrt(1/8) sqrt(1/8) -1/4 1/4 0 -sqrt(1/32) sqrt(9/8) -sqrt(1/32)
sqrt(1/8) sqrt(1/8) -1/4 3/4 0 0 -sqrt(1/8) -sqrt(1/8)
sqrt(1/8) sqrt(1/8) -1/4 3/4 0 0 -sqrt(1/8) sqrt(9/8)
];

Summary of programs:

1. outputcdf53.m – main function output CDF 5/3 transform matrix
2. legall53td.m – forward CDF 5/3 wavelet transform
3. ilegall53td.m – inverse CDF 5/3 wavelet transform

download zip

Leave a Reply

Your email address will not be published. Required fields are marked *