www.gusucode.com > wavelet 源码程序 matlab案例代码 > wavelet/LiftingExample.m
%% Lifting % This example shows how to use lifting on a 1-D signal. %% % Create a 1-D signal that is piecewise constant over 2 samples. Add $N(0,0.1^{2})$ % noise to the signal. x = [1 1 2 2 -3.5 -3.5 4.3 4.3 6 6 -4.5 -4.5 2.2 2.2 -1.5 -1.5]; x = repmat(x,1,64); rng default; x = x+ 0.1*randn(size(x)); %% % Plot the signal and zoom in on the first 100 samples to visualize the % correlation in neighboring samples. plot(x); set(gca,'xlim',[0 100]); %% % Use the _lazy_ wavelet to obtain the even and odd polyphase components % of the signal. LS = liftwave('lazy'); [A,D] = lwt(x,LS); %% % If you plot the detail (wavelet) coefficients in |D|, you see that this % transform has not decorrelated the signal. The wavelet coefficients look % very much like the signal. %% % Add a dual lifting step that subtracts the even-indexed coefficient from % the odd-coefficient one sample later, $x(2n+1)-x(2n)$. els = {'d',-1,0}; LSnew = addlift(LS,els); %% % Because the signal is piecewise constant over consecutive samples with % additive noise, the new dual lifting step should result in wavelet coefficients % small in absolute value. In this case, the wavelet transform does decorrelate % the data. Verify this by finding the approximation and detail coefficients % with the new dual lifting step. [A,D] = lwt(x,LSnew); %% % If you plot the detail (wavelet) coefficients, you see that the wavelet % coefficients no longer resemble the original signal. %% % The approximation coefficients, |A|, of the previous transform constitute % the even polyphase component of the signal. Therefore, the coefficients % are affected by aliasing. Use a primal lifting step to update the approximation % coefficients and reduce aliasing. The primal step replaces the approximation % coefficients by $x(2n)+1/2(x(2n+1)-x(2n))$, which is equal to the average % of $x(2n)$ and $x(2n+1)$. The averaging is a lowpass filtering, which % helps to reduce aliasing. els = {'p',1/2, 0}; LSnew = addlift(LSnew,els); %% % Use the updated lifting scheme to obtain the wavelet transform of the % input signal. [A,D] = lwt(x,LSnew); %% % Add the appropriate scaling to ensure perfect reconstruction. Obtain the % approximation and wavelet coefficients using lifting and reconstruct the % signal using |ilwt|. Verify perfect reconstruction. LSnew(end,:) = {sqrt(2),sqrt(2)/2,[]}; [A,D] = lwt(x,LSnew); x1 = ilwt(A,D,LSnew); max(abs(x1-x)) %% % The preceding example designed a wavelet, which effectively removed a % zero-th order polynomial (constant). If the behavior of the signal is % better represented by a higher-order polynomial, you can design a dual % wavelet with the appropriate number of vanishing moments to decorrelate % the signal. %% % Use the lifting scheme to design a wavelet with 2 vanishing moments. A % dual wavelet with 2 vanishing moments decorrelates a signal with local % behavior approximated by a first-order polynomial. Create a signal characterized % by first-order polynomial behavior with additive $N(0,0.25^{2})$ noise. y = [1 0 0 4 0 0 -1 0 0 2 0 0 7 0 0 -4 0 0 1 0 0 -3]; x1 = 1:(21/1024):22-(21/1024); y1 = interp1(1:22,y,x1,'linear'); rng default; y1 = y1+0.25*randn(size(y1)); plot(x1,y1); set(gca,'xlim',[1 22]); %% % In this case, the wavelet coefficients should remove a first-order polynomial. % If the signal value at an odd index, $x(2n+1)$, is well approximated by % a first-order polynomial fitted to the surrounding sample values, then % $1/2(x(2n)+x(2n+2))$ should provide a good fit for $x(2n+1)$. In other % words, $x(2n+1)$ should be the midpoint between $x(2n)$ and $x(2n+2)$. % % It follows that $x(2n+1)-1/2(x(2n)+x(2n+2))$ should decorrelate the signal. % % Start with the lazy wavelet transform and add a dual lifting step which % models the preceding equation. LS = liftwave('lazy'); els = {'d',[-1/2 -1/2],1}; LSnew = addlift(LS,els); %% % Use the lifting scheme to obtain the approximation and detail coefficients % and plot the result. [A,D] = lwt(y1,LSnew); subplot(211) plot(A); set(gca,'xlim',[1 512]); title('Approximation Coefficients'); subplot(212) plot(D); set(gca,'xlim',[1 512]); title('Detail Coefficients'); %% % You see that the wavelet coefficients appear to only contain noise, while % the approximation coefficients represent a denoised version of the original % signal. Because the preceding transform uses only the even polyphase component % for the approximation coefficients, you can reduce aliasing by adding % a primal lifting step. Finally, add the normalization constants to produce % a perfect reconstruction filter bank. %% % Obtain the discrete wavelet transform with the new lifting scheme and % plot the results. els = {'p',[1/4 1/4],0}; LSnew = addlift(LSnew,els); LSnew(end,:) = {sqrt(2),sqrt(2)/2,[]}; [A,D] = lwt(y1,LSnew); subplot(211) plot(A); set(gca,'xlim',[1 512]); title('Approximation Coefficients'); subplot(212) plot(D); set(gca,'xlim',[1 512]); title('Detail Coefficients'); %% % Demonstrate that you have designed a perfect reconstruction filter bank. y2 = ilwt(A,D,LSnew); max(abs(y2-y1))