小波分析實(shí)驗(yàn)報(bào)告

班級(jí)：電信0906小波分析實(shí)驗(yàn)?zāi)夸浶〔ǚ治鰧?shí)驗(yàn)一2.二維圖像的傅里葉頻譜分析3.一維采樣信號(hào)的傅里葉頻譜分析4.小波的分解與重構(gòu)(a)隨機(jī)離散信號(hào)的分解與重構(gòu)(b)二維圖像的分解與完全重構(gòu)小波分析實(shí)驗(yàn)二(a)低頻保留壓縮(c)全局閾值壓縮(d)基于Birge-Massart策略中經(jīng)驗(yàn)系數(shù)的分層閾值壓縮(e)基于小波包變換的圖像壓縮(f)基于JPEG2K標(biāo)準(zhǔn)的圖像壓縮1.2JPEG與JPEG2K的綜合比較(a)小波基的選擇(b)小波變換系數(shù)分析2.基于小波的盲數(shù)字水印2.2DCT與DWT盲水印綜合比較小波分析實(shí)驗(yàn)一(1)實(shí)驗(yàn)圖像Figure1.(1)原始信號(hào)及其Fourier頻譜(連續(xù)表示)Figure1.(2)原始信號(hào)及其Fourier頻譜(離散表示)(2)實(shí)驗(yàn)代碼clears=rand(1,NUM)*MAXN;subplot(2,1,1),plot(s);xlabel('discreten');ylabel('magnituden');subplot(2,1,2),Sf=abs(fft(s));xlabel('discreten');ylabel('magnitude');figure,subplot(2,1,1),stem(s,'*');xlabel('discreten');ylabel('magnituden');subplot(2,1,2),stem(Sf,'*');MAXY=max(Sf)xlabel('discreten');ylabel('magnitude');title('Fourier頻譜(離散表示));代碼實(shí)現(xiàn)對(duì)長(zhǎng)度為32的隨機(jī)一維行矩陣進(jìn)行傅里葉變換，并顯2、對(duì)一幅離散二維圖像(如512*512像素尺度的黑色背景疊加一個(gè)20*40像素尺寸的白色矩形，也可選擇其他圖像，如fingerprint等，可使用imread函數(shù)讀取圖像),實(shí)現(xiàn)傅里葉變換，并顯示原圖(1)實(shí)驗(yàn)圖像Figure2.(1)原始圖像及其Fourier頻譜(一)Figure2.(2)原始圖像及其Fourier頻譜(2)實(shí)驗(yàn)代碼號(hào)對(duì)一幅離散二維圖像實(shí)現(xiàn)傅里葉變換，并顯示原圖及傅里葉頻譜S=fftshift(log(1+abs(F));S=gscale(S);f=zeros(512,512);f(247:266,237:276)=1;S=abs(F):S=gscale(S);S=abs(Fc);S=gscale(S);subplot(2,2,3),imshow(S);title('居中的傅里葉頻譜圖像');S=fftshift(log(1+abs(F)));S=gscale(S);(3)實(shí)驗(yàn)分析實(shí)驗(yàn)圖像中的第二個(gè)子圖是對(duì)原始圖像直接進(jìn)行傅里葉變換的結(jié)果，圖像四角上的亮點(diǎn)是圖像周期性導(dǎo)致的結(jié)果。為了使變換的原點(diǎn)移到頻率矩形的中心，我是用了函數(shù)fftshift,新產(chǎn)生的圖像如第三個(gè)子圖所示，可以看出居中的效果還是很明顯的。不過(guò)，雖然該移動(dòng)nagnitudephasenagnitudephase如期望的那樣完成了，但頻譜中值的動(dòng)態(tài)范圍與8比特顯示(此時(shí)中心處的明亮值占支配地位)相比要大得多。于是我用對(duì)數(shù)變換對(duì)原始圖像進(jìn)行處理，效果如第四個(gè)子圖所示，可見(jiàn)可視細(xì)節(jié)的增加還是很明顯的。另外在圖像的現(xiàn)實(shí)中，我還調(diào)用了自編的m函數(shù)gscale來(lái)實(shí)現(xiàn)按比例調(diào)整輸入圖像的亮度，以求達(dá)到最佳的可視化效果。函數(shù)gscale是被標(biāo)度的圖像，[0,255],后者將輸出圖像標(biāo)度為全范圍[0,65535]。3、對(duì)一個(gè)連續(xù)的函數(shù)(如8cos(18x)+2cos(x)等)進(jìn)行采樣，得到離散信號(hào)，然后實(shí)現(xiàn)傅里葉變換，顯示原信號(hào)及傅里葉頻譜。(1)實(shí)驗(yàn)圖像原始信號(hào)3+8*原始信號(hào)3+8*cos(18*t+pi/6)+2*cos(t+pi/4)傅立葉頻譜(連續(xù)表示)5000葉actualfrequency相位(弧度)頻率曲線(xiàn)圖5actualfrequency400800900600名原始信號(hào)3+8*cos(18*t+pi/6)+2*cos(t+pi/4)clear號(hào)確定采樣區(qū)間0...Tx=zeros(1,N);y=zeros(1,N);x(k)=(k/N)*T;y(k)=3+8*cos(18*x(k)+pi/6)+2*cos(x(k)+pi/4)figure,subplot(4,1,1),plot(x,y);title('原始信號(hào)3+8*cos(18*t+pi/6)+2*cos(t+pi/4)');xlabel('discreten');ylabel('magnitude');subplot(4,1,2),xlabel('discreten');ylabel('magnitude');subplot(4,1,3),xlabel('actualfrequency');ylabel('magnitude');title('幅度—頻率曲線(xiàn)圖');8P(i)=P(i)*180/pi;g轉(zhuǎn)換為角度(3)實(shí)驗(yàn)分析代碼實(shí)現(xiàn)將原始信號(hào)3+8*作為實(shí)驗(yàn)的拓展部分，我還繪制了幅度—頻率曲線(xiàn)圖和相位(弧度)—頻率曲線(xiàn)圖。在幅度—頻率曲線(xiàn)圖中將離散頻率點(diǎn)用公式angle函數(shù)計(jì)算相位，得到頻率成分對(duì)應(yīng)的相位，并用弧度表示(如果函數(shù)),顯示原信號(hào)、分解后的小波子帶、重構(gòu)(1)實(shí)驗(yàn)圖像Figure4.(1)Haar小波一維重構(gòu)與分解Figure4.(2)Haar小波一維重構(gòu)與分解各層系數(shù)ReconstructedsignalFigure4.(6)Db2小波一維重構(gòu)與分解各層系數(shù)(2)實(shí)驗(yàn)代碼x=rand(1,16)*MAXN;[lpr,hpr]=wfilters('haar','r');[lpd,hpd]=wfilters('haar','d');wavdanr1d(x,lpd,hpd,lpr,hpr,dim);[lpr,hpr]=wfilters('db2','r');[lpd,hpd]=wfilters('db2','d');wavdanr1d(x,Ipd,hpd,Ipr,hpr,dim)[lpr,hpr]=wfilters('db3','r');[lpd,hpd]=wfilters('db3','d');wavdanr1d(x,Ipd,hpd,lpr,hpr,dim)實(shí)驗(yàn)結(jié)果顯示：原始信號(hào)與重構(gòu)信號(hào)的均方根誤差為：原始信號(hào)與重構(gòu)信號(hào)的均方根誤差為：原始信號(hào)與重構(gòu)信號(hào)的均方根誤差為：(3)實(shí)驗(yàn)分析函數(shù)分別實(shí)現(xiàn)小波的一級(jí)分解和一級(jí)重構(gòu)，其中waveletdec為失真完全重構(gòu)，這也是小波變換最大的特點(diǎn)之一。自帶的小波分解重構(gòu)函數(shù)DWT2,wavedec2函數(shù)),顯示原圖、分解后的小波子帶、重構(gòu)后的(1)實(shí)驗(yàn)圖像Figure4.(8)Haar小波分解與重構(gòu)圖像db3重構(gòu)圖像clearimg=imread('house.png');img=im2double(img);figure(1),[lpdhpd]=wfilters('haar','d');[lprhpr]=wfilters('haar','r');imgr=wavrec2(Il,hl,Ih,hh,Ir,Ic,Ipr,hpr);figure(2),subplot(1,2,1);subplot(1,2,2);rmse1=compare(img,imgr);disp('Haar小波重構(gòu)圖像均方根誤差為：');lpdhpd]=wfilters('db2'.'d');[lprhpr]=wfilters('db2','r');imgr=wavrec2(ll,hl,Ih,hh,Ir,Ic,Ipr,hpr);figure(3),subplot(1,2,1);subplot(1,2,2);rmse2=compare(img,imgr);disp('Db2小波重構(gòu)圖像均方根誤差為：');[lprhpr]=wfilters('db3','r');rmse3實(shí)驗(yàn)結(jié)果顯示：rmse1=1.6726e-016rmse2=5.4105e-014rmse3=8.2068e-013(3)實(shí)驗(yàn)分析數(shù)。其中upspl,downspl,waveletdec和waveletrec函數(shù)功能與小波對(duì)一幅圖像實(shí)現(xiàn)分解(在本實(shí)驗(yàn)中，僅進(jìn)行1級(jí)分解),而wavdec2函數(shù)則在waveletrec基礎(chǔ)上用指定小波對(duì)一幅變換后的圖像實(shí)現(xiàn)重全重構(gòu)，這也是小波變換最大的特點(diǎn)之一。作為該實(shí)驗(yàn)的更高級(jí)拓展，我運(yùn)用Gonzales提供的關(guān)于小波快個(gè)子函數(shù)和一個(gè)主函數(shù)。主函數(shù)為experiment4b2函數(shù)。子函數(shù)是wavefilter,wavefast,wavework,wavecut,wavecopy,wavepaste,waveback,Figure4.(12)Haar小波圖像分解(細(xì)節(jié)系數(shù)放大8倍Figure4.(13)Haar小波圖像分解(細(xì)節(jié)系數(shù)絕對(duì)值放大8倍)Figure4.(14)Haar小波重構(gòu)圖像Figure4.(16)Db2小波圖像分解(細(xì)節(jié)系數(shù)放大8倍)db2分解(細(xì)節(jié)系數(shù)絕對(duì)值放大8倍)Figure4.(17)Db2小波圖像分解(細(xì)節(jié)系數(shù)絕對(duì)值放大8倍Figure4.(18)Db2小波重構(gòu)圖像Figure4.(20)Db3小波圖像分解(細(xì)節(jié)系數(shù)放大8倍)Figure4.(21)Db3小波圖像分解(細(xì)節(jié)系數(shù)絕對(duì)值放大8倍)Figure4.(22)Db3小波重構(gòu)圖像實(shí)驗(yàn)代碼%對(duì)一幅灰度圖像，實(shí)現(xiàn)用Haar小波、Db2、Db3小波的分解與f=imread('lena.png');f=imread('lena.png');f=imread('lena.png');(c)對(duì)一幅灰度圖像(可使用house,lena圖像等)用Haar小波或Db小波進(jìn)行分解，對(duì)小波系數(shù)做修改(如閾值等),最后再重構(gòu)。將原(1)實(shí)驗(yàn)圖像重構(gòu)圖像2clearimg=imread('house.png');img=im2double(img);figure(1),[lpdhpd]=wfilters('db4','d');[lprhpr]=wfilters('db4','r');alpha=0.25;maxcoef=max(imgd(:));mincoef=min(imgd(:));vth=alpha*(maxcoef-mincoef);hl(find(Il<vth))=0;Ih(find(Il<vth))=0;hh(find(Il<vth))=0;imgr=wavrec2(II,hl,Ih,hh,Ir,Ic,Ipr,hpr);figure(2),subplot(1,2,1);imshow(imgd);title('分解圖像1);subplot(1,2,2);rmse1=compare(img,imgr);disp('Db4小波重構(gòu)圖像均方根誤差為(閾值一):');maxcoef=max(imgd(:));mincoef=min(imgd(:));vth=alpha*(maxcoef-mincoef);rmse1=0.0256rmse2=0.0356(3)實(shí)驗(yàn)分析在該題中，我在experiment4_b基礎(chǔ)上，對(duì)分解的高頻系數(shù)進(jìn)行對(duì)高頻部分，重構(gòu)的圖像與原始圖像在視覺(jué)感知上可我們也可以利用這一點(diǎn)進(jìn)行圖像的壓縮(關(guān)于基于小波的圖像壓縮，F(xiàn)igure4.(28)第一級(jí)細(xì)節(jié)系數(shù)置零重構(gòu)實(shí)驗(yàn)代碼%對(duì)一幅灰度圖像(可使用house,lena系數(shù)做修改(如閾值等),最后再重構(gòu)圖像等)用Haar小波或Db小波進(jìn)行分解，對(duì)小波clearf=imread('house.png');figure,[c1,s1]=wavefast(f,4,'figure,在這個(gè)實(shí)驗(yàn)中，我利用函數(shù)wavefast得到小波分解向量c和記錄矩陣s,緊接著依次將第一層、第二層、第三層、第四層小波分解的高頻系數(shù)置零，得到的圖像分別如因?yàn)樗菍⒃瓐D中第一層細(xì)節(jié)系數(shù)置零而得到的。figure4.(29)的模糊程度更明顯了，因?yàn)樗鼘⒌诙蛹?xì)節(jié)系數(shù)置零，而第二層細(xì)節(jié)系數(shù)中劇的模糊程度。小波分析實(shí)驗(yàn)實(shí)驗(yàn)四基于小波的圖像壓縮請(qǐng)參考“小波分析及其應(yīng)用”課件用簡(jiǎn)單方案實(shí)現(xiàn)基于小波變換的圖像壓縮。算法一：只保留低頻部分，舍去高頻圖像的重要信息和絕大部分能量都集中在低頻。而大部分高頻系數(shù)都接近于0,而且頻率越高這種現(xiàn)象越明顯。于是去掉高頻部分LH,HL,HH,就成為了一種最簡(jiǎn)單有效地方案，還可以對(duì)LL進(jìn)行多次分解，去掉更多的高頻成分。(1)實(shí)驗(yàn)圖像Figure4.(1)原始圖像Figure4.(2)Bior3.7小波分解圖像Figure4.(3)第一次壓縮圖像Figure4.(4)第二次壓縮圖像(2)實(shí)驗(yàn)代碼title('原始圖像);disp('壓縮前圖像X的大小為：);c1=[a1,h1;v1,d1];8下面進(jìn)行圖像壓縮處理%保留圖像分解第一層低頻信息，進(jìn)行圖像壓縮ca1=appcoef2(c,s,'bior3.7',1);ca1=wcodemat(ca1,440,'mat',0);ca1=0.5*ca1;figure(3),image(ca1);colormap(map);title('第一次壓縮');disp('第一次壓縮圖像的大小為：');whos('ca1')g保留圖像分解第二層低頻信息，進(jìn)行圖像壓縮號(hào)第二層的低頻信息即為ca2,顯示第二層低頻信息g首先對(duì)第二層進(jìn)行量化編碼ca2=appcoef2(c,s,'bior3.7',2);ca2=wcodemat(ca2,440,'mat',0);ca2=0.25*ca2;figure(4),image(ca2);colormap(map);title(‘第二次壓縮');disp('第二次壓縮圖像大小為：');whos('ca2')壓縮前圖像X的大小為：NameSizeBytesAttributes第一次壓縮圖像的大小為NameBytesAttributes第二次壓縮圖像的大小為：NameBytes45000Attributes縮圖像可以看出隨著分解級(jí)數(shù)的增加，壓縮的圖像逐漸變得模糊，出現(xiàn)較大損失。所以在算法一中，根據(jù)預(yù)期的壓縮效果和適當(dāng)?shù)囊曈X(jué)效果而確定可靠的分解級(jí)數(shù)是最理應(yīng)事先考慮的。由于小波分析固有的時(shí)頻特性，我們可以在時(shí)頻兩個(gè)方向?qū)π〔ㄏ禂?shù)進(jìn)行處理，這樣就可以對(duì)我們感興趣的部分提供不同的壓縮精于是我們可以把圖像中部的細(xì)節(jié)系數(shù)置零，而圖像其余部分的細(xì)節(jié)系數(shù)以及所有的低頻系數(shù)均保持不變，已達(dá)到預(yù)先期待的壓縮部分細(xì)節(jié)而保留部分細(xì)節(jié)的目的。(1)實(shí)驗(yàn)圖像Figure4.(5)原始圖像Figure4.(6)壓縮圖像Figure4.(7)一層分解后個(gè)系數(shù)圖像Figure4.(6)處理后各層系數(shù)圖像(2)實(shí)驗(yàn)代碼使用sym4小波對(duì)信號(hào)進(jìn)行一層小波分解合復(fù)制原圖像小波系數(shù)rca1=ca1;rch1=ch1;%將3個(gè)細(xì)節(jié)系數(shù)的中間部分置零rch1(33:97,33:97)=zeros(65,65);rcv1(33:97,33:97)=zeros(65,65);codrca1=wcodemat(rca1,192);codrch1=wcodemat(rch1,192);codrcv1=wcodemat(rcv1,192);S將處理后的系數(shù)圖像組合為一個(gè)圖像codrx=[codrca1,codrch1;codrcv1,codrcd1]。充值處理后的數(shù)據(jù)figure(1);image(wcodemat(X,192));colormap(map);8求壓縮信號(hào)的能量成分disp('壓縮信號(hào)的能量成分比為：');disp('壓縮信號(hào)與原信號(hào)的標(biāo)準(zhǔn)差為：');壓縮信號(hào)的能量成分為：per=1.0000(3)實(shí)驗(yàn)分析算法二用于實(shí)現(xiàn)對(duì)感興趣區(qū)域和不感興趣區(qū)域?qū)崿F(xiàn)不同精度的可以看出只有中間部分變得模糊(如在原圖中很清晰的圍巾的條紋不能分辨),而其他部分的細(xì)節(jié)信息仍然很清楚。(1)實(shí)驗(yàn)圖像Figure4.(9)原始圖像(2)實(shí)驗(yàn)代碼Figure4.(10)壓縮圖像%顯示圖像實(shí)驗(yàn)結(jié)果顯示：perf12=99.9784(3)實(shí)驗(yàn)分析由實(shí)驗(yàn)圖像可知，全局壓縮在一定的壓縮程度內(nèi)，可以得到可視化效果較好的壓縮圖像。算法四：基于Birge-Massart策略中經(jīng)驗(yàn)系數(shù)的分層閾值壓縮(1)實(shí)驗(yàn)圖像能量成分999858%零系數(shù)成分49.8893%Figure4.(12)全局閾值壓縮圖像Figure4.(13)分層閾值壓縮壓縮圖像(2)實(shí)驗(yàn)代碼求得顏色映射表的長(zhǎng)度，以便后面的轉(zhuǎn)換用默認(rèn)的方式求出圖像的全局閾值8顯示全局閾值系數(shù)號(hào)對(duì)圖像做全局閾值用bior.3.5小波對(duì)圖像進(jìn)行3層分解指定Birge-Massart策略中的經(jīng)驗(yàn)系數(shù)alpha=1.5;m=2.7*prod(s(1,:));。根據(jù)各層小波系數(shù)確定分層閾值8對(duì)原圖像做分層閾值[xd1,cxd1,sxd1,perf01,perfl21]=wdencmp('Ivd',c,s,'bior3.5',3,thr1,'s')thr1figure(1);figure(2);image(wcodemat(xd,nbc));colormap(map);disp('全局壓縮中置零的系數(shù)個(gè)數(shù)百分比';figure(3);image(wcodemat(xd1,nbc));colormap(map);disp('分層閾值壓縮中置零的系數(shù)個(gè)數(shù)百分比');perfl21實(shí)驗(yàn)結(jié)果顯示：全局閾值系數(shù)為：thr=4.0000thr1=47.343347.343347.343366.188166.188166.1881全局壓縮中置0的系數(shù)個(gè)數(shù)百分比：perf0=49.8693全局壓縮后圖像剩余能量百分比：分層閾值壓縮中置0的系數(shù)個(gè)數(shù)百分比：(3)實(shí)驗(yàn)分析該實(shí)驗(yàn)中，我同時(shí)使用了求默認(rèn)閾值的ddencmp命令和基于經(jīng)驗(yàn)基于小波包變換的圖像壓縮(1)實(shí)驗(yàn)圖像Figure4.(14)原始圖像Figure4.(15)小波包分解全局閾值化壓縮圖像(2P2123P1242R#RTER風(fēng)86RkPapsFigure4.(16)最有小波樹(shù)(2)實(shí)驗(yàn)代碼。求顏色索引表長(zhǎng)度nbc=size(map,1);。得到閾值系數(shù)，保留層數(shù)和小波樹(shù)優(yōu)化標(biāo)準(zhǔn)。通過(guò)以上得到的參數(shù)對(duì)圖像進(jìn)行壓縮。更改索引表為pink索引表壓縮后圖像剩余能量百分比：(3)實(shí)驗(yàn)分析R,R,R,+1R,+1算法六：用jpeg2000標(biāo)準(zhǔn)實(shí)現(xiàn)圖像壓縮用jpeg2000對(duì)圖像進(jìn)行壓縮處理的算法：首先將原圖像的灰度值減去2m-1來(lái)進(jìn)行圖像的灰度級(jí)移動(dòng)。然后可計(jì)算圖像的每一行和每一列的小波變換，對(duì)圖像的每一層的子帶進(jìn)行量化操作，子帶b的系數(shù)ab(u,v)可利用公式得到量化值：量化步長(zhǎng)△b對(duì)于無(wú)損壓縮為1,對(duì)于有損壓縮沒(méi)有規(guī)定具體的計(jì)算規(guī)則，但是可以按照下面的方法進(jìn)行：同，若RI是圖像的位深則可以按照下圖方式計(jì)算y,為能量權(quán)(energyweight),表示量化時(shí)變換系數(shù)產(chǎn)生的誤差以及由此而引起的圖像誤差之間的關(guān)系，也就是變換系數(shù)的單位誤差的平方差的和。編碼處理的最后步驟是在比特平面的基礎(chǔ)上對(duì)量化后的系數(shù)進(jìn)行算術(shù)編碼。(1)實(shí)驗(yàn)圖像Figure4.(5)原始圖像Figure4.(6)JPEG2000Figure4.(7)差異圖像Figure4.(8)均方根誤差直方圖(2)實(shí)驗(yàn)代碼8用JPEG2000實(shí)現(xiàn)圖像壓縮compressed_value=imratio(f,c1);計(jì)算圖像壓縮比f(wàn)igure,figure,實(shí)驗(yàn)結(jié)果顯示：compressed_value圖像均方根誤差為：rms=3.9792(3)實(shí)驗(yàn)分析在這里涉及到兩個(gè)核心函數(shù)分別是：im2jpeg2k和jpeg2k2im。前者實(shí)現(xiàn)的功能是將圖像進(jìn)行壓縮編碼得到一個(gè)結(jié)構(gòu)體變量c1這里面含有圖像編碼的一些信息。而jpeg2k2im函數(shù)實(shí)現(xiàn)的功能是解碼將顯示重現(xiàn)圖像。JPEG與JPEG2K的綜合比較：(1)壓縮圖像比較Figure4.(9)原始圖像Figure4.(12)JPEG差異圖像Figure4.(13)JPEFigure4.(14)JPEG均方根誤差直方圖Figure4.(15)JPEG2000均方根誤差直方圖(2)壓縮結(jié)果比較compressed_value=15.8664compressed_value=37.3318=2.6251rms=3.9792(3)綜合比較分析為成熟的工業(yè)化標(biāo)準(zhǔn)，壓縮效果和壓縮后圖像的可視化效果較算法是真正基于需求而進(jìn)行的利用多種其他技術(shù)(特別是數(shù)據(jù)的編碼和解碼算法)的圖像壓縮算法。由此可見(jiàn)，利用小波分析進(jìn)行的圖像壓縮基于小波圖像壓縮技術(shù)的算法研究(一)小波基的選擇(1)小波基的正則性和消失矩(2)小波基的線(xiàn)性相位(3)所處理圖像與小波基的相似性(4)小波函數(shù)的能量集中性(5)綜合考慮壓縮效率和計(jì)算復(fù)雜度程度。正則性對(duì)圖像壓縮有一定的影響，如果圖像大部分是光滑的，(1)實(shí)驗(yàn)圖像Figure4.(16)原始圖像尺度為1時(shí)的低領(lǐng)圖像Figure4.(17)Haar第一層壓縮圖像尺度為2時(shí)的低頻圖像2002Figure4.(18)sym8第一層壓縮圖像尺度為2時(shí)的低新圖像50100150200250Figure4.(19)Haar第二層壓縮圖像Figure4.(20)sym8第二層壓縮圖像Figure4.(21)Haar第三層壓縮圖像Figure4.(22)sym8第三層壓縮圖像(2)實(shí)驗(yàn)代碼gray=imread('cameraman.tif');subplot(2,2,1),image(wcodemat(A1,nbcol));colormap(map);title('尺度為1時(shí)的低頻圖像');subplot(2,2,2),image(wcodemat(H1,nbcol));colormap(map);subplot(2,2,3),image(wcodemat(V1*90,nbcol));colormap(map);subplot(2,2,4),image(wcodemat(D1*89,nbcol));colormap(map);figure(2)subplot(2,2,1),image(wcodemat(A2,nbcol));colormap(map);title('尺度為1時(shí)的低頻圖像');subplot(2,2,2),image(wcodemat(H2,nbcol));colormap(map);subplot(2,2,3),image(wcodemat(V2,nbcol));colormap(map);subplot(2,2,4),image(wcodemat(D2,nbcol));colormap(map);figure(3)subplot(2,2,1),image(wcodemat(A3,nbcol));colormap(map);title('尺度為1時(shí)的低頻圖像');subplot(2,2,2),image(wcodemat(H3,nbcol));colormap(map);subplot(2,2,3),image(wcodemat(V3,nbcol));colormap(map);subplot(2,2,4),image(wcodemat(D3,nbcol));colormap(map);gray=imread('cameraman.tif');[X,map]=gray2ind(gray);figure(4);figure(1)subplot(2,2,1),image(wcodemat(A1,nbcol));colormap(map);title('尺度為1時(shí)的低頻圖像');subplot(2,2,2),image(wcodemat(H1,nbcol));colormap(map);subplot(2,2,3),image(wcodemat(V1*90,nbcol));colormap(map);subplot(2,2,4),image(wcodemat(D1*89,nbcol));colormap(map);figure(2)subplot(2,2,1),image(wcodemat(A2,nbcol));colormap(map);title('尺度為1時(shí)的低頻圖像');subplot(2,2,2),image(wcodemat(H2,nbcol));colormap(map);subplot(2,2,3),image(wcodemat(V2,nbcol));colormap(map);subplot(2,2,4),image(wcodemat(D2,nbcol));colormap(map);figure(3)subplot(2,2,1),image(wcodemat(A3,nbcol));colormap(map);title('尺度為1時(shí)的低頻圖像');image(wcodemat(D3,nbcol));colormap(map);(3)實(shí)驗(yàn)分析(二)小波變換系數(shù)分析密度函數(shù)(PDF),小波變換后各子代概率密度函數(shù)可以通過(guò)統(tǒng)計(jì)方法逼近。資料表明，在高頻子代，小波系數(shù)更符合廣義高斯分布。下面的實(shí)驗(yàn)對(duì)小波系數(shù)分布理論進(jìn)行了驗(yàn)證。(1)實(shí)驗(yàn)圖像Figure4.(23)尺度為1時(shí)的系數(shù)分布Figure4.(24)尺度為2時(shí)的系數(shù)分布Figure4.(25)尺度為3時(shí)的系數(shù)分布(2)實(shí)驗(yàn)代碼8%小波分解window=[77];Opname='mean2';IsIH1=IH1./(sigmlH1+.1^10);IsIV1=IV1./(sigmIV1+.1^10);IsID1=ID1./(sigmlD1+.1^10);IsIH2=IH2./(sigmlH2+.1^10);IsIV2=IV2./(sigmIV2+.1^10);IsID2=ID2./(sigmID2+.1^10);IsIH3=IH3./(sigmIH3+.1^10);IsIV3=IV3./(sigmIV3+.1^10);IsID3=ID3./(sigmID3+.1^10);Iss1=reshape(Isl11,S1(1)*S1(2),1);Iss2=reshape(IsIH1,S1(1)*S1(2),1);Iss3=reshape(IsIV1,S1(1)*S1(2),1);Iss4=reshape(IsID1,S1(1)*S1(2),1);S2=size(lsl21);Iss21=reshape(Isl21,S2(1)*S2(2),1);Iss22=reshape(IsIH2,S2(1)*S2(2),1);Iss23=reshape(IsIV2,S2(1)*S2(2),1);Iss24=reshape(IsID2,S2(1)*S2(2),1);S3=size(Isl31);title((b)尺度為1時(shí)的水平高頻系數(shù)分布);subplot(2,2,3),histfit(Iss3);title((c)尺度為1時(shí)的垂直高頻系數(shù)分布);尺度為2時(shí)的水平高頻系數(shù)分布);histfit(lss31);尺度為3時(shí)的低頻系數(shù)分布);尺度為3時(shí)的水平高頻系數(shù)分布);histfit(lss33);(3)實(shí)驗(yàn)分析波對(duì)原圖進(jìn)行3次分解后得到的低頻和高頻圖像。Figure4.(23)、figure4.(24)、figure4.(25)是對(duì)figure4.(18)、figure4.(20)、figure4.(22)所對(duì)應(yīng)的小波系數(shù)的分差應(yīng)該服從正態(tài)分布，小波系數(shù)的均值為零。Figure4.(23)、可以看出，在尺度為1的情況下，低通小波系數(shù)不服從正態(tài)分布而高實(shí)驗(yàn)五基于小波的盲數(shù)字水印(1)實(shí)驗(yàn)圖像watermarkwatermarkFigure5.(1)原始圖像Figure5.(2)水印圖像Figure5.(3)嵌入水印圖像Figure5.(4)提取水印圖像高斯低通濾波高斯低通濾波旋轉(zhuǎn)90度部分剪切旋轉(zhuǎn)90度部分剪切Figure5.(5)水印圖像進(jìn)行攻擊直接提取水印加入高斯白噪聲提取水印高斯低通濾波提取水印5|JPEG壓縮提取水印部分剪切提取水印旋轉(zhuǎn)90度提取水印Figure5.(6)攻擊后的圖像提取水印(2)實(shí)驗(yàn)代碼%攻擊實(shí)驗(yàn)測(cè)試魯棒性image1=image1+noise;M2=image2;image3=imread('watermarkedM3=image3;image4=imwtwm;image4(1:128,1:128)=256;M4=image4;image5=imwtwm;image5=imrotate(image5,90,'crop');image5=mat2gray(image5);M5=image5;M6=imwtwmsubplot(2,3,1),subplot(2,3,2),subplot(2,3,3),subplot(2,3,4),subplot(2,3,5),subplot(2,3,6),figure(6),wmext1=wmextract(imwtwm,0.2,2);wmext2=wmextract(image1,0.2,2);wmext3=wmextract(image2,0.2,2);wmext4=wmextract(image3,0.2,2);wmext5=wmextract(image4,0.2,2);wmext6=wmextract(image5,0.2,2);subplot(2,3,1),imshow(wmext1),title('直接提取水印);subplot(2,3,2),imshow(wmext2),title('加入高斯白噪聲提取水印')subplot(2,3,3),(3)實(shí)驗(yàn)分析1.用Daubechies雙正交9/7小波對(duì)圖像進(jìn)行3級(jí)分解。得到低頻分量的小波系數(shù)f(LL3,x,y)、水平分量小波系數(shù)f(HLn,x,y)和對(duì)角分量小波系數(shù)f(HHn,x,y),n=1,2,32.文字水印信息用置亂算法進(jìn)行置亂(對(duì)水印圖像進(jìn)行加密)3.將置亂后的文字水印w(x,y)一級(jí)小波分解：低頻分量小波系數(shù)w(LL,x,y)水平分量小波系數(shù)w(LH,x,y)垂直分量小波系數(shù)w(HL,x,y)對(duì)角分量小波系數(shù)w(HH,x,y)4.對(duì)小波變化后的原始圖像的第二級(jí)的高頻子帶和第三級(jí)高頻F’(x,y)=F(x,y)+aw(x,y)原始圖像(x,y)的位置上嵌入的水印小波系數(shù)值1.對(duì)嵌入水印圖像f’(x,y)進(jìn)行三級(jí)小波分解，得到低頻分量小波系數(shù)f’(HLn,x,y)對(duì)角分量小波系數(shù)f’(HHn,x,y),n=1,2,3同理對(duì)原始圖像f(x,y)進(jìn)行3級(jí)小波分解2.對(duì)小波變化后的含水印圖像和原始圖像進(jìn)行提取算法處理：w(x,y)=(F’(x,y)-F(x,y))/a波系數(shù)4.將上一步得到的水印信息進(jìn)行逆置亂變換(即解密),得到提取作為實(shí)驗(yàn)的拓展部分，我用5種方法對(duì)嵌入水印后的圖像進(jìn)行了NontensorProductWaveletFilterBanks[J].IEEETransactionsonimageprocessing,VOL.19,NO.12,2010.3271-3284MISIFigure5.(7)原始圖像Figure5.(8)水印圖像旋轉(zhuǎn)90度旋轉(zhuǎn)90度Figure5.(10)水印圖像進(jìn)行攻擊加入高斯白噪聲提取水印高斯低通濾波提取水印加入高斯白噪聲提取水印高斯低通濾波提取水印L部分剪切提取水印部分剪切提取水印clearl=imread('watermarkjpg);T=graythresh(I);l=im2bw(I,T);figure(1),diml=size(I);rm=diml(1);cm=diml(2);alpha=50;k1=randn(1,8)k2=randn(1,8);aO=imread('lena.jpg');aO=rgb2gray(a0);pnsr_cover=double(a0);figure(2),[r,c]=size(a0);cdaO=blkproc(a0,[8,8],'dct2');cda1=cda0;fori=1:rmforj=1:cmx=(i-1)*8;y=(j-1)*8;elsek=k2:endcda1(x+1,y+8)=cda0(x+1,y+8)+alpha*k(1)cda1(x+2,y+7)=cda0(x+2,y+7)+alpha*k(2);cda1(x+3,y+6)=cda0(x+3,y+6)+alpha*k(3);cda1(x+4,y+5)=cda0(x+4,y+5)+alpha*k(4);cda1(x+5,y+4)=cda0(x+5,y+4)+alpha*k(5);cda1(x+6,y+3)=cda0(x+6,y+3)+alpha*k(6);cda1(x+7,y+2)=cda0(x+7,y+2)+alpha*k(7);cda1(x+8,y+1)=cda0(x+8,y+1)+alpha*k(8);%5旋轉(zhuǎn)90度%start_time=cuptime;image1=ima