Open in the MATLAB editor: c09fd_example
function c09fd_example fprintf('c09fd example results\n\n'); m = int64(7); n = int64(6); fr = int64(5); wavnam = 'Bior1.1'; mode = 'period'; wtrans = 'Multilevel'; a = zeros(m, n, fr); a(:, :, 1) = [3, 2, 2, 2, 1, 1; 2, 9, 1, 2, 1, 3; 2, 5, 1, 2, 1, 1; 1, 6, 2, 2, 7, 2; 5, 3, 2, 2, 4, 7; 2, 2, 1, 1, 2, 1; 6, 2, 1, 3, 6, 9]; a(:, :, 2) = [2, 1, 5, 1, 2, 3; 2, 9, 5, 2, 1, 2; 2, 3, 2, 7, 1, 1; 2, 1, 1, 2, 3, 1; 2, 1, 2, 8, 3, 3; 1, 4, 5, 1, 2, 7; 8, 1, 3, 9, 1, 2]; a(:, :, 3) = [3, 1, 4, 1, 1, 1; 1, 1, 2, 1, 2, 6; 4, 1, 7, 2, 5, 6; 3, 2, 1, 5, 9, 5; 1, 1, 2, 2, 2, 1; 2, 6, 3, 9, 5, 1; 1, 1, 8, 2, 1, 3]; a(:, :, 4) = [5, 8, 1, 2, 2, 1; 1, 2, 2, 9, 2, 9; 2, 2, 2, 1, 1, 3; 1, 1, 1, 5, 1, 2; 3, 2, 8, 1, 9, 2; 2, 1, 9, 1, 2, 2; 3, 6, 5, 3, 2, 2]; a(:, :, 5) = [5, 2, 1, 2, 1, 1; 3, 1, 9, 1, 2, 1; 2, 3, 1, 1, 7, 2; 7, 2, 2, 6, 1, 1; 5, 1, 7, 2, 1, 1; 2, 1, 3, 2, 2, 1; 5, 3, 9, 1, 4, 1]; % Query wavelet filter dimensions [nwl, nf, nwct, nwcn, nwcfr, icomm, ifail] = ... c09ac(wavnam, wtrans, mode, m, n, fr); % Perform Discrete Wavelet transform [c, dwtlvm, dwtlvn, dwtlvfr, icomm, ifail] = c09fc(n, fr, a, nwct, nwl, icomm); fprintf(' Number of Levels : %d\n\n', nwl); fprintf(' Number of coefficients in 1st dimension for each level:\n'); fprintf(' %8d', dwtlvm(1:nwl)); fprintf('\n'); fprintf(' Number of coefficients in 2nd dimension for each level:\n'); fprintf(' %8d', dwtlvn(1:nwl)); fprintf('\n'); fprintf(' Number of coefficients in 3rd dimension for each level:\n'); fprintf(' %8d', dwtlvfr(1:nwl)); fprintf('\n'); % Print the first level HLL coefficients want_level = 1; % Select the approximation coefficients. want_coeffs = 4; % Identify each set of coefficients in c for ilevel = nwl:-1:1 if ilevel ~= want_level continue end nwcm = dwtlvm(nwl-ilevel+1); nwcn = dwtlvn(nwl-ilevel+1); nwcfr = dwtlvfr(nwl-ilevel+1); fprintf('\n--------------------------------\n'); fprintf(' Level %d output is %d by %d by %d.\n', ilevel, nwcm, nwcn, nwcfr); fprintf('--------------------------------\n\n'); for itype_coeffs = 0:7 if itype_coeffs ~= want_coeffs continue end % Unless we're looking at the deepest level of nesting, which contains % approximation coefficients, advance the pointer on past the preceding % levels if ilevel == nwl locc = 0; else locc = 8*dwtlvm(1)*dwtlvn(1)*dwtlvfr(1); for i = ilevel + 1 : nwl - 1 locc = locc + 7*dwtlvm(nwl-i+1)*dwtlvn(nwl-i+1)*dwtlvfr(nwl-i+1); end end % Now decide which coefficient type we are considering switch (itype_coeffs) case {0} if (ilevel==nwl) fprintf('Approximation coefficients (LLL)\n'); locc = locc + 1; end case {1} fprintf('Detail coefficients (LLH)\n'); if (ilevel==nwl) % Advance pointer past approximation coefficients locc = locc + nwcm*nwcn*nwcfr + 1; else locc = locc + 1; end case {2} fprintf('Detail coefficients (LHL)\n'); if (ilevel==nwl) % Advance pointer past approximation coefficients and 1 set of % detail coefficients locc = locc + 2*nwcm*nwcn*nwcfr + 1; else % Advance pointer past 1 set of detail coefficients locc = locc + nwcm*nwcn*nwcfr + 1; end case {3} fprintf('Detail coefficients (LHH)\n'); if (ilevel==nwl) % Advance pointer past approximation coefficients and 2 sets of % detail coefficients locc = locc + 3*nwcm*nwcn*nwcfr + 1; else % Advance pointer past 2 sets of detail coefficients locc = locc + 2*nwcm*nwcn*nwcfr + 1; end case {4} fprintf('Detail coefficients (HLL)\n'); if (ilevel==nwl) % Advance pointer past approximation coefficients and 3 sets of % detail coefficients locc = locc + 4*nwcm*nwcn*nwcfr + 1; else % Advance pointer past 3 sets of detail coefficients locc = locc + 3*nwcm*nwcn*nwcfr + 1; end case {5} fprintf('Detail coefficients (HLH)\n'); if (ilevel==nwl) % Advance pointer past approximation coefficients and 4 sets of % detail coefficients locc = locc + 5*nwcm*nwcn*nwcfr + 1; else % Advance pointer past 4 sets of detail coefficients locc = locc + 4*nwcm*nwcn*nwcfr + 1; end case {6} fprintf('Detail coefficients (HHL)\n'); if (ilevel==nwl) % Advance pointer past approximation coefficients and 5 sets of % detail coefficients locc = locc + 6*nwcm*nwcn*nwcfr + 1; else % Advance pointer past 4 sets of detail coefficients locc = locc + 5*nwcm*nwcn*nwcfr + 1; end case {7} fprintf('Detail coefficients (HHH)\n'); if (ilevel==nwl) % Advance pointer past approximation coefficients and 6 sets of % detail coefficients locc = locc + 7*nwcm*nwcn*nwcfr + 1; else % Advance pointer past 5 sets of detail coefficients locc = locc + 6*nwcm*nwcn*nwcfr + 1; end end if itype_coeffs > 0 || ilevel == nwl if (itype_coeffs==0) % For a multi level transform approx coeffs stored as % nwcm x nwcn x nwcfr i1 = locc; for k = 1:nwcfr for j = 1:nwcn for i = 1:nwcm d(i,j,k) = c(i1); i1 = i1 + 1; end end end else % ... but detail coefficients are stored as ncwfr x nwcm x nwcn for k = 1:nwcfr for j = 1:nwcn for i = 1:nwcm i1 = locc - 1 + (j-1)*nwcfr*nwcm + (i-1)*nwcfr + k; d(i,j,k) = c(i1); end end end end % Print out the selected set of coefficients fprintf('Level %d, Coefficients %d:\n', ilevel, itype_coeffs); for k = 1:nwcfr fprintf('Frame %d:\n', k); for i = 1:nwcm for j=1:nwcn fprintf('%8.4f ', d(i, j, k)); end fprintf('\n'); end end end end end % Reconstruct original data [b, ifail] = c09fd(nwl, c, m, n, fr, icomm); % Check reconstruction matches original eps = 10*double(m*n*fr)*x02aj; err = a-b; frob = 0; for i=1:fr fnew = sqrt(sum(sum(err(:,:,i).^2))); frob = max(frob,fnew); end if frob < eps fprintf('\nSuccess: the reconstruction matches the original.\n'); else fprintf('\nFail: Frobenius norm of b-a is too large.\n'); end
c09fd example results Number of Levels : 2 Number of coefficients in 1st dimension for each level: 2 4 Number of coefficients in 2nd dimension for each level: 2 3 Number of coefficients in 3rd dimension for each level: 2 3 -------------------------------- Level 1 output is 4 by 3 by 3. -------------------------------- Detail coefficients (HLL) Level 1, Coefficients 4: Frame 1: -4.9497 0.0000 0.0000 0.7071 1.7678 -3.1820 0.7071 2.1213 1.7678 0.0000 0.0000 0.0000 Frame 2: 4.2426 -2.1213 -4.9497 0.7071 -0.0000 -0.7071 -1.4142 -3.1820 1.4142 0.0000 0.0000 0.0000 Frame 3: 2.1213 -4.9497 -0.7071 -2.8284 -4.2426 4.9497 2.1213 2.8284 -0.7071 0.0000 0.0000 0.0000 Success: the reconstruction matches the original.