小波实验报告双树复小波变换

《小波实验报告双树复小波变换》由会员分享，可在线阅读，更多相关《小波实验报告双树复小波变换（19页珍藏版）》请在装配图网上搜索。

1、一.题目：双树复小波变换二、目的：双树复小波和实小波变换的比较 算法及其实现：提取阶梯型边界点1. 算法。= 0” (0 + J 匕、 3),#) =几必+ J(j,久) 幅值：I 叽卫)|= Jd心,) +/(/,)相位：Zd (j, n) = arctailUhUJ2. 代码实现。我采用 Matlab 函数编程实现。具体程序见 shift_test_2D.m, diawciic.m setfig.m, dnvavexfin2in, dtwaveifin2.m,waveifiii2.m.wavexfiii2.mSkelMap.nio咚 0,爲和屮hWh分别是双正交对偶尺度函数与对偶小波，力0

2、（），人）（），力1（丹） 斤（/?）是相应的低通滤波器和高通滤波器，即它们满足实部：0力二工九（）0力-“）虎部：0&=血工goOMg川）屮k （0 =忑工入（“加-）匕=忑工& 0必 -刃）n必（f） =工 hQ（77） 一 ）冈(/) = 工 /?, (/?M (2r - n)MZDE-町II讥（0 =迈工ii 0诫 -双树复小波变换可以通过离散小波变换DWT实现：一个DWT产生实部，另一个 产生虚部。四、实现工具：Matlab五、程序代码：(1) shift_test_2D.nl：% shift_test_2D.m% M-flle to peiiorm a 4level wavelet

3、 transform on a circle using Q-shift % dual wavelet tree and DWT, and to compaie sliift iiivaiiance properties. % Nick Kingsbury、Cambridge University. May 2002 clear allclose all% Draw a ciicular disc.x = round(drawcirc(64 J,0,0,256) - 0.5) * 200); setfig(l);colorniap(gray(256) iniage(niHi(niax(x+12

4、8,1),256);set(gca/position0.1 0.25 .25 .5);axis(,off);axis(,miage,);% draw(xx);tiqinput (256 x 256)*FontSize14); diawnow% Do 4 levels of CWTYLYh = dtwavexfiii2(x,4/neai_sym_b,qsluft_b,);% Loop to feconstiuct output from coefs at each level ill turn.% Starts with tlie finest level.titl = flst,;,2ndf;

5、3rdt;,4th,;,Low,J;vy = zeros(size(x)* 2 3); vtl = l:size(x,l); yt2 = l:size(x,2);for nilev = 1:5,mask = zeros(6,5);mask(:anlev) = 1;z = dtwaveifiii2(Yl*mask( 1,5), Yh/neai_sym_bVqshift_bniask); figure; draw(z) ;drawnowyy(yti,yt2)= z；yt2 = vt2 + size(x,2)/2;end% dispfPress a key% pause% Now do same w

6、ith DWT.% Do 4 levels of Real DWT using fantomni, (9,7)-tap filters. YLYh = wavexfiii2(x,4；antonuir);vtl = l:size(x,l) + size(x.l); yt2 = l:size(x,2);for nilev = 1:5,mask = zeros(3,5);mask(:anlev) = 1;z = waveifiii2( Yl*niask( 1,5). Yli/antoninimask);figure; draw(z) ;drawnowyy(ytl,yt2) = z;yt2 = vt2

7、 + size(x,2)/2;endfiguie;setfig(gcf);colorniap(gray(256)iniage(niHi(niax(yv+128,1)256);set(gca,position0.1 0.1 .8 .8);axisfoff);axis(,iinage,);hold onplot(128*l;l*l:4 0;6+l,128*0;4*l 1 1 1 2；2+l；-k*)；hold offdle(Coinponents of reconstmcted disc” iniages*FontSize 14); text(-0.01*size(yv7).025*size(yy

8、4)?DT CWTThonzJ); text(0.02*size(yy,2)4.02*size(yyj),vravelets：7horiz,r7vert,/t,);text(-0.01 *size(yy,2),0.75 *size(yyJ),DWT7horiz7rl)；for k=l:4, text(k* 128-63,size(yy, 1)*1.02,sprintf(level %dk)/FontSize14/honz,/c,/vertVf); endtext(5 * 128+1 ,size(yy, 1)*1.02,level 4 scalmg fii.,FontSize14/horiz7c

9、lvertVt,);diawnow% print -deps cucrecq.epsdisp(,Piess a key to see perfect reconstruction propeity .)pause% Accumulate the images from lowband upwards to show perfect reconstmction. sy = size(x,2)/2;fbrnilev = 4l:l,yt2 = l:sy + (nilev-l)*sy;yy(:,yt2) = yy(:,yt2) + yy(:,yt2+sy);endfigure；setfig(gcf);

10、colorniap(gray(256)iniage(niHi(niax(yv+128,1)256);set(gca,position0.1 0.1 .8 .8);axisfoff);axis(,iinage,);title(,Accumulated leconstnictions fiom each level of DT CWT TontSize14); text(size(yy,2)*0.5,size(vyJ)*1.02/Accumulated reconstmctions fiom each level of DWT TontSize 14；hor7cVertyt1);diawnowre

11、turn(2) functionp = diawcirc(i;wJx.dy,N)% function p = diawcirc(r,w.dx,dy；N)% Generate an image of size N*N pels, contauung a circle% radius r pels and centred at dx.dy relative% to tlie centre of tlie image The edge of the ciicle is a cosine shaped% edge of width w (from 10 to 90% pomts).x = ones(N

12、J) * (1:N(N+l)/2 - dx)/r);y = (kN* - (N+l)/2 - dv)/i) * ones(l,N);p = 0.5 + 0.5 * sin(iniii(inax(exp(-0.5 * (x +y) exp(O.5)*(T*3/w), -pi/2), pi/2);return(3) functionZ = dnvaveifiii2(YLYh.biort,qshift,gain_mask);% Function to peiibrm an n-level dual-tree complex wavelet (DTCWT)% 2-D reconstmction.% Z

13、 = dnvaveifiii2(YLYh,bioit,qshift,gain_mask);%Y1 - The real lowpass image from the filial level%Yh - A cell array contauiiiig the 6 complex highpass subiinages for each level.% bion - ntonini* = Aiitonini 9,7 tap filters%legair = LeGall 5,3 tap filters.%,neai_sym_a, = Near-Synunetric 5,7 tap filters

14、.%,neai_sym_b, = Near-Synunetric 13,19 tap filters% qshift - qshifL06= Quarter Sample Sluft Onhogonal (Q-Sluft) 10.10 tap filters,%(only 6,6 non-zero taps)%qshift_cf = Q-sliift 10,10 tap filters,%(with 10,10 non-zero taps, unlike qshift_06)%qshift_b =Q-Sliift14,14 tap filters.%qshift_c=Q-Sliift16,16

15、 tap filters.%qshift_d =Q-Sliift18,18 tap filters.% gam_mask - Gain to be applied to each subband.%gain_mask(d、l) is gain for subband with diiection d at level 1.%If gain_mask(d)=0, no computation is performed for band (d)%Default gain_mask = ones(6Jength(Yh).% Z - Reconstmcted real image matrix% Fo

16、i example: Z = dmraveifiii2(YLYh?neai_sym_b,qshift_b,);% performs a 3level reconstmction fiom Yl.Yh using the 13.19-tap filters% for level 1 and the Q-shift 14-tap filters for levels = 2.% Nick Knigsbury and Ciaii Shaffiey% Cainbiidge Univeisity, May 2002a = length(Yh); % No of levels.if nargin = 2;

17、 %tliis ensures that for level -1 we never do the following111 = c2q(YhcuiTentJevel(:/,! 6),gaiii_mask(l 6,current_level);hl = c2q(YhcuiTent_level(:，:,3 4),gaiii_mask(3 4,current_level);lili = c2q(Ylicurrent_level(:/,2 5),gam_mask(2 5,cunent_level);% Do even Qshift filters on columns, yl = colifilt(

18、 Z,g0b,g0a) + colifilt(lh,g 1 b,g 1 a);y2 = colifilt(hl,gOb,gOa) + colifilt(hh,g 1 b,g 1 a);% Do even Qshift filters on lows.Z = (colifilt(y 1 /,gOb,gOa) + colifilt(y2. g 1 b,g 1 a).% Check size of Z and ciop as requuediow_size coLsize = size(Z);S = 2*size(Yhcunent_level-l);if iow_size = S(l) %check

19、 to see if tliis result needs to be cropped for the rows Z = Z(2:iow_size-1,:);endif coLsize = S(2) %check to see if tliis result needs to be cropped for the cols Z = Z(:,2:col_size-1);endif any(size(Z) = S(l:2)，enoiSizes of subbands are not valid for DTXVtWEIFM2r);endcunentjevel = current_level -1;

20、 end if cunentjevel = 1;111 = c2q(YhcuiTentJevel(:,:,! 6),gaiii_mask(l 6,cuiTent_level); hl = c2q(YhcuiTentJevel(：,：,3 4),gaiii_mask(3 4,cuiTent_level); lili = c2q(Ylicurrent_level(:,:,2 5),gam_mask(2 5,cunent_level);% Do odd top-level filters on colunms.yl = colfilter( Z,gOo) + colfilter(lh,glo);y2

21、 = colfilter(lil,gOo) + colfilter(hh,glo);% Do odd top-level filters on rows.Z = (colfilter(yl.gOo) + colfilter(y2 .g 1end return%=%* INTERNAL FUNCTION *%= fiinction x = c2q(w,gam)% function z = c2q(w,gain)% Scale by gain and convert fiom complexto real quad-numbers in 乙% Arrange pixels from tlie re

22、al and irnag pans of the 2 subbands% into 4 separate subiniages % A一一BReIinofw(:,:,l)% II% II% C一一DReImofw(j，2)sw = size(w);x = zeros(2*sw(l:2);if anv(w(:) & any(gaiii)sc = sqrt(O5) * gain;P = w(:,:,l)*sc(l) + w(:,:,2)*sc(2);Q = w(:,:,l)*sc -w(:,:,2)*sc(2);tl = l:2:size(x,l);t2 = l:2:size(x,2);% Rec

23、over each of the 4 cornei s of the quadsx(tl,t2) = real(P); % a = (A+C)*sc;x(tl,t2+l)= iinag(P); % b = (B+D)*sc;x(tl+l,t2)= iinag(Q); % c = (B-D)*sc;x(tl+l,t2+l) = -real(Q); % d = (C-A)*sc;endreturn(3) function YLYh.Yscale = dtvavexfiii2(X,nlevels,biort,qsluft);% Function to peiibrm a n-level DTCWT-

24、2D decompostion on a 2D matrix X% YLYhYscale = dhvavexfiii2(X,iilevels,bion.qshift);% X - 2D real matrix-linage% iilevels - No. of levels of wavelet decomposition% bion - antoninr = Aiitonini 9,7 tap filters%legall1 = LeGall 5,3 tap filters.%,neai_sym_a, = Near-Symmetric 5,7 tap filters.,neai_sym_b,

25、 = Near-Synunetric 13,19 tap filters% qshift - qshifL06= Quarter Sample Sluft Onhogonal (Q-Sluft) 10.10 tap filters, %(only 6.6 noil-zero taps)%qshifta* = Q-sliift 10,10 tap filters,%(with 10,10 non-zero taps, unlike qshift_06)%qshift_b = Q-Sliift14,14 tap filters.%qshift* = Q-Sliift16,16 tap filter

26、s.%qshiftd* = Q-Sliift1&18 tap filters.%Y1- The real lowpass image from the filial level%Yh- A cell array contammg the 6 complex highpass subiniages for each level.% Yscale - Tliis is an OPTIONAL output aigument, that is a cell array contaniing%real lowpass coefficients for eveiy scale% % Example: Y

27、LYh = dnvavexfiii2(X,3,neai_sym_b;% performs a 3level transform on the real image X using the 13.19-tap filters% for level 1 and the Q-shift 14-tap filters for levels = 2.% % Nick Kuigsbuiy and Ciaii Shaffiey % Cainbiidge University, Sept 2001if isstr(biort) & issti(qshifl) %Check if the mputs are s

28、trings biort_exist = exist(biort .mat);qshift_exist = exist(qshift r.maf);if biort_exist = 2 & qshift_exist = 2;%Check to see if the inputs existas niat filesload (bioit);load (qshift);elseenoiXPlease enter the correct names of the Bioithogonal or Q-Sluft Filters, see help DTWAVEXFM2 for details/);e

29、ndelseeuoPlease enter tlie names of the Bioithogonal or Q-Shift Filters as shown ill help DTWAVEXFM2J;endorginaLsize = size(X);if nduns(X) = 3;euospimtThe entered image is %dx%dx%d. please enter each image slicesepaiatelyorginal_size(l),orgmal_size(2),orgiiial_size(3);end % The next few lines of cod

30、e check to see if tlie image is odd m size, if so an extra % low/colunm will be added to the bottonVright of the imageinitial row extend = 0; %iiiitialiseinitial col extend = 0;if aiiy(rem(orgiiiaLsize( 1 ),2), %if sx(l) is not divisable by 2 then we need to extend X by adding a row at tlie bottomX

31、= X; X(end,:);% Any fiirther extension will be done in due course.initial row extend = 1:endif any(rem(orgmaLsize(2),2), %if sx(2) is not divisable by 2 then we need to extend X by addmg a col to the leftX= X X(:,end);%Any fbrther extension will be done in due courseinitial col extend = 1;endextende

32、d_size = size(X);if iilevels = 0. return; end%initialiseYh=cell(iilevels, 1);if nargout = 3Yscale=cell(nlevels, 1);%tliis is only requiied if tlie user specifies a thiid outputcomponent.endS = ;sx = size(X);if iilevels = 1,% Do odd top-level filters on cols.Lo = colfilter(X,hOo).r;Hi = colfilter(X,h

33、lo).r;% Do odd top-level filters on rows.LoLo = colfiltei(Lo.hOo).f;% LoLoYh 1 = zeios(size(LoLo)/2 6);% Horizontal pair% Vertical pair% Diagonal paifYhl(:二1 6) = q2c(colfilteT(Hhh0o)；Yhl(:二3 4) = q2c(colfilteT(Loilo)J;Yhl(:二2 5) = q2c(colfiltei(HLhlo)/);S = size(LoLo) ;S;if naigout = 3Yscalel = LoL

34、o;if iilevels = 2;for level = 2:nlevels;row_size col_size = size(LoLo);if any(rem(row_size,4). % Extend by 2 rows if no. of rows of LoLo aie divisableLoLo = LoLo(l,:); LoLo: LoLo(end,:);endif any(rem(coLsize,4).% Extend by 2 cols if no. of cols of LoLo are divisable by4；LoLo = LoLo(:,l) LoLo LoLo(:,

35、end);end % Do even Qshift filters on rowsLo = coldfilLoLoJiObJiOa).1;Hi = coldfilt(LoLo,h 1 b,h 1 a).1;% Do even Qshift filters on columns.LoLo = coldfilt(LoJiOb、hOa)： %LoLoYhlevel = zeros(size(LoLo)/2 6);Yhlevel(:,:,! 6) = q2c(coldfilt(Hi,h0b,h0a) J; % Horizontal Yhlevel(:,:,3 4) = q2c(coldfilt(Lo,

36、hlbJ】la) J; % Vertical Yhlevel(:,:,2 5) = q2c(coldfilt(Hi,hlb,h 1 a).1);% DiagonalS = size(LoLo) ;S;if naigout = 3Yscale level = LoLo;endendendY1 = LoLo;if iiiitiaLrow_extend = 1 & initial_col_extend = 1;wainmg(sprintf(, 1V The image entered is now a %dx%d NOT a %dx%d r The bottom row and rightmost

37、column have been duplicated, prior to decomposition. 1V ；extended_size( l),extended_size(2),orgmaLsize( l),oigiiial_size(2);end if mitial row extend = 1 ;wainmg(sprintf(, rr The image entered is now a %dx%d NOT a %dx%d r Row number %d has been duplicated, and added to the bottom of the linage, piior

38、 to decomposition. extended_size( 1 ),extended_size(2),orgiiial_size( 1 ),orginaLsize(2),orginal_size(l);end if initial col extend = 1;warnmg(spiintf(, rr The miage entered is now a %dx%d NOT a %dx%d r Col number %d has been duplicated, and added to the right of the image, prior to decomposition. ex

39、tended_size( l),extended_size(2),orgiiial_size( l),orginaLsize(2),orginal_size(2);endreturn%* INTERNAL FUNCTION * fiinction z = q2c(y) % function z = q2c(v)% Convert from quads in y to complex numbers in 乙sv = size(y)；tl = l:2:sy(l)； t2 = l:2:sy(2);j2 = sqrt(0.5 -0.5);% Arrange pixels from tlie corn

40、ers of the quads into % 2 subiinages of alternate real and irnag pixels % ab% II% II% c-d% Combine (久b) and (d,c) to form two complex subiinages. p = y(tl,t2)*j2(l) + y(tl,t2+l)*j2(2);% p = (a + jb) / sqrt(2)q = y(tl+l,t2+l)*j2(l) - y(tl+l,t2)*j2(2); % q = (d -jc)/ sqrt(2) % Form the 2 subbands in 乙

41、z = cat(3.p-q,p+q);letuin(4) functionZ = dnvaveifiii2(YLYh.biort,qshift,gain_mask);% Function to peiiorm an n-level dual-tree complex wavelet (DTCWT)% 2-D reconstmction.% Z = dtwaveifiii2(YLYh,bioit,qshift,gain_mask);%Y1 - The real lowpass image from the filial level%Yh - A cell array containing the

42、 6 complex highpass subiinages for each level.% bion - antonini* = Aiitonini 9,7 tap filters.%legair = LeGall 5,3 tap filters.%,neai_sym_a, = Near-Symmetric 5,7 tap filters.%,neai_sym_b, = Near-Synmietric 13,19 tap filters% qshift - qshift_06= Quarter Sample Sluft Onhogonal (Q-Sluft) 10.10 tap filte

43、rs, %(only 6.6 noil-zero taps)%qshifta* = Q-sliift 10,10 tap filters,%(with 10,10 non-zero taps, unlike qshift_06)%qshift_b = Q-Sliift 14,14 tap filters.%qshift_c= Q-Sliift 16,16 tap filters.%qshiftd1 = Q-Sliift 1&18 tap filters.% gain_mask - Gain to be applied to each subband.%gain_inask(d、l) is ga

44、in for subband with diiection d at level 1.%If gam_mask(dj) = 0. no computation is performed for band (dj).%Default gain_mask = ones(6Jength(Yh).% Z - ReconstRicted real image matrix% Foi example: Z = dmraveifiii2(YLYh/neai_syni_b,qshift_b,);% performs a 3level reconstmction fiom Yl.Yh using the 13,

45、19-tap filters% for level 1 and the Q-shift 14-tap filters for levels = 2.% Nick Kuigsbury and Ciaii Shaffiey % Cainbiidge University, May 2002 a = length(Yh); % No of levels.if nargin = 2; %tliis ensures that for level -1 we never do the following111 = c2q(YhcuiTentJevel(:,:,! 6),gaiii_mask(l 6,cur

46、rent_level);hl = c2q(YhcuiTentJevel(：,：,3 4),gaiii_mask(3 4,current_level);lili = c2q(Ylicurrent_level(:,:,2 5),gam_mask(2 5,cunent_level);% Do even Qshift filters on columns.yl = colifilt( Z，gOb，gOa) + colifilt(lh,g 1 b,g 1 a);y2 = colifilt(hl,gOb,gOa) + colifilt(hli,g 1 b,g 1 a);% Do even Qshift f

47、ilters on lows.Z = (colifilt(yl .gOb,gOa) + coiifilt(y2.,5glb,gla).,j% Check size of Z and crop as requuedrow_size col_size = size(Z);S = 2*size(Yhcunent_level-l);if iow_size = S(l) %check to see if tliis result needs to be cropped for the rows Z = Z(2:iow_size-1,:);endif col_size = S(2) %check to s

48、ee if tliis result needs to be cropped for the cols Z = Z(:,2:col_size-1);endif any(size(Z) = S(l:2),eiTor(fSizes of subbands are not valid for DTXVtWEIFM2f);endcunenOevel = cunentjevel -1; endif cunenOevel = 1;111 = c2q(YhcuirenClevel(:,:,! 6),gaiii_mask(l 6,current_level); hl = c2q(YhcuirenClevel(

49、：,：,3 4),gaiii_mask(3 4,current_level); lili = c2q(Ylicurrent_level(:,:,2 5),gam_mask(2 5,cunent_level);% Do odd top-level filters on colunms. yl = colfilter( Z,gOo) + colfilter(lh,glo);y2 = colfilter(lil,gOo) + colfilter(hli,glo);% Do odd top-level filters on rows.Z = (colfilter(yl.gOo) + colfilter

50、(y2 .g 1 o).*;end leturn%=%* INTERNAL FUNCTION *%= fiinction x = c2q(w.gam)% function z = c2q(w,gain)% Scale by gain and convert fiom complex w(j,l:2) to real quad-numbers in 乙% Arrange pixels from tlie real and miag pans of the 2 subbands% into 4 separate subiniages % A一一BReIinofw(:,:,l)% II% II% C

51、一一DReImofw(j，2)sw = size(w);x = zeros(2*sw(l:2);if anv(w(:) & any(gaiii)sc = sqrt(O5) * gain;P = w(:,:,l)*sc(l) +w(:,:,2)*sc ;Q = w(:,:,l)*sc(l) - w(:,:,2)*sc(2);tl = l:2:size(x,l);t2 = l:2:size(x,2);% Recover each of the 4 corners of the quadsx(tl,t2) = real(P); % a = (A+C)*sc;x(tl,t2+l)= iinag(P);

52、 % b = (B+D)*sc;x(tl+l,t2)= iinag(Q); % c = (B-D)*sc;x(tl+l,t2+l) = -real(Q); % d = (C-A)*sc;return(5) fiinction YLYhYscale = wavexfiii2(Xjilevels,biort);% Function to peiibrm a n-level DWT-2D decompostion on a 2-D matrix X.% % YLYhYscale = dhvavexfiii2(Xailevels,biort);% % X - real 1-D signal colum

53、n vector (or matrix of vectors)% iilevels - No. of levels of wavelet decomposition% % biort - antonini* = Aiitonini 9,7 tap filters%legalF = LeGall 5,3 tap filters.%,neai_sym_a, = Near-Synmietric 5,7 tap filters.%,neai_sym_b, = Near-Synmietric 13,19 tap filters.%Y1- The lowpass subband from tlie fin

54、al level.%Yh- A cell array contammg the highpass subband for each level.% Yscale - Tins is an OPTIONAL output aigument, that is a cell airay contammg%the lowpass coefficients at eveiy scale% % Example: YLYli = wavexfiii2(X.4/near_sym_b,);% performs a 4level 2D DWT on the real image X using the 13,19

55、-tap filters.%if isstr(biort) % Check if tlie bioit input is a stringbiort_exist = exist(biort .mat);if biort_exist = 2, % Check to see if the filter exists as a mat file load (bioit);elseenoiXPlease enter the correct name of the Bioithogonal Filter, see help WAVEXFM2 for details/);endelseenor(*Plea

56、se enter the name of the Biorthogonal Filter as shown in help WAVEXFM2.);endL = size(X);if anv(rem(L,2), % ensuie tliat X is an even length, thus enabling it to be extended if needs be. error(*Size of X must be a multiple of 21);end%initialiseYh=cell(iilevels, 1);if nargout = 3Yscale=cell(nlevels, 1);% Tliis is only requiied if the user specifies a thud outputcomponent.endLoLo = X;for level = 1 :iilevels;if iem(size(LoLo, 1)4),% Check to see if height of LoLo is divisable by 4, if not extend.LoLo = LoLo(l,:); LoLo; LoLo(end/);endif