圖像數字水印+matlab程序

1、第三章圖像數字水印的方案3.1 圖像數字水印的技術方案在數據庫中存儲在國際互聯網上傳輸的水印圖像一般會被壓縮，有時達到很高的壓縮比。因此，數字水印算法所面臨的第一個考驗就是壓縮。JPEG和EZW（EmbeddedZero-TreeWavelet）壓縮是最常見的兩種壓縮方法。JPEG是基于離散余弦變換域的壓縮方法，而EZW是基于小波變換域的壓縮方法。前人的研究證明采用與壓縮算法相同的變換域水印方法，對于壓縮的穩(wěn)健性較強。因此，我研究圖像文件水印算法主要集中在變換域算法及利用人眼視覺特性上。數字水印的嵌入要求即要考慮視覺透明性，又要保證嵌入水印后圖像的穩(wěn)健性，這兩個方面存在著矛盾。保證視覺透明性，

2、就要將水印嵌入到人眼不敏感區(qū)，也就是嵌入到圖像的高頻分量中。而多數圖像處理方法對于圖像高頻部分的損壞程度較高，如有損壓縮、高頻濾波等。水印很容易在經歷圖像處理的過程中丟失。這樣，則無法保證圖像數字水印的穩(wěn)健性。如果要獲得很好的穩(wěn)健性，數字水印應加在人眼敏感的低頻部分，圖像的大部分能量集中在低頻部分，如果對于低頻部分進行處理，水印固然會失去，而圖像也沒有了利用價值，然而，水印的嵌入會對圖像的質量有非常大的影響，這又無法保證視覺透明性。數字水印算法的實現基本分為三個部分：宿主圖像的變換，水印的嵌入和水印的檢測，分別描述如下。3.2 基于DCTM的圖像數字水印技術離散余弦變換（DiscreteCos

3、ineTransform）屬于正交變換圖像編碼方法中的一種。正交變換圖像編碼始于1968年。當時安德魯斯（Andrews）等人發(fā)現大多數自然圖像的高頻分量相對幅度較低，可完全舍棄或者只用少數碼字編碼，提出不對圖像本身編碼，只對其二維傅立葉（DFT1系數進行編碼和傳輸。但DFT是一種正交變換，運算量很大，常常使實時處理發(fā)生困難，第二年他們就用Walsh-Hadamard變換（WHT）取代DFT可以使運算量明顯減少，這是因為WHT變換只有加減法而無需乘法。但是更有意義的是離散余弦變換和離散正旋變換的出現，它們具有快速算法，精確度高。其中最重要的是1974年提出的DCT,因為其變換矩陣的基向量很近似

4、于托伯利茲矩陣的特征向量，而托伯利茲矩陣又體現了人類語言及圖像信號的相關性。因此，DCT常常被認為是語音與圖像信號變換的準最佳變換。圖像是二維的，所以在研究時主要用到二維DCT,以及二維IDCT來對圖像進行處理。3.2.1 離散余弦變換（DCT的定義數字圖像X（m,n）是具有M行N列的一個矩陣。為了同時減弱或去除圖像數據相關性，可以運用二維DCT,將圖像從空間域轉換到DCT變換域。根據定義，二維離散余弦變換（DCT）定義如下：2M1NY(k,l)：c(k)c(l)%X(m,n)cos,：/MNm2n乂(2m1)k二2Mcos(2nF二2Nc(l)式中：m,k=0,1,M1;n,l=0,1,N。

5、其中函數c(k)=Jq2k=01k=1,2,，M-1J1/行k=0=1k=1,2,，N-1二維逆離散余弦變換(IDCT)的定義如下：2MN(2m1)k二(2n1)l二X(m,n)=、c(k)cY(k,l)coscos(.MNKw2M2N式中：m,k=0,1,M1;n,l=0,1,N1。3.2.2 離散余弦變換的特點在基于DCT的變換編碼中，圖像是先經分塊（8X8或16X16）后再經DCT,這種變換是局部的，只反映了圖像某一部分的信息。當然也可以對整幅圖像的特點，但是運算速度比分塊DCT要慢。圖像經DCT后，得到的DCT圖像有三個特點：一是系數值全部集中到0值附近（從直方圖統計的意義上），動態(tài)范

6、圍很小，這說明用較小的量化比特數即可表示DCT系數；二是DCT變換后圖像能量集中在圖像的低頻部分，即DCT圖像中不為零的系數大部分集中在一起（左上角），因此編碼效率很高。三是沒有保留原圖像塊的精細結構，從中反映不了原圖像塊的邊緣、輪廓等信息，這一特點是由DCT缺乏時局域性造成的。如下左圖31是原始圖彳t經過DCT變換后的系數圖像為圖32。兩條線劃分出圖像的低頻、中頻和高頻分別所在的矩形區(qū)域。可以看出，圖像DCT變換后大部分參數接近于零，只有左上角的低頻部分有較大的數值，中頻部分參數值相對較小,而大部分圖31:原圖像圖32:變換后的系數圖像3.2.3 離散余弦變換的數字水印算法根據離散余弦變換后

7、的參數性質，本文采用了以ZigZag方式重排變換域系數的方法，選出中頻分量，用數字水印序列對其進行非線性調制。水印檢測時，待檢測圖像仍按比方式選擇變換域系數，與待水印進行相關運算，與閾值比較來判斷是否所含水印。離散余弦域的數字水印算法的具體實現分為三步：宿主圖像的變換，數字水印的嵌入，數字水印的檢測。3.2.3.1 宿主圖像的DC儂換對于NXN大小的256灰度級的宿主圖像I進彳TNXN二維離散余弦變換(DCT)。以ZigZag方式對于DCT變換后的圖像頻率系數重新排列成一維向量Y=yi,y2,nxn.并取出序列中第L+1到L+M的中頻系數部分彳導到YL=Yl+i,YL+2,lYm3.2.3.2

8、 數字水印的嵌入假設數字水印w為一服從標準正態(tài)分布的隨機實數序列，用數字序列表示為W=W1,W2,Wmo用W對Y序列中第L+1到L+M的中頻系數部分的值進行修改，按以下公式進行：y=Vi控制系統(controlsystem)神經網絡(neuralnetwork)、圖像處理(imageprocessing)魯棒控制(robustcontrol)、非線性系統控制設計(nonlinearsystemcontroldesign)、系統辨識(systemidentification)、最優(yōu)化(optimization)、模糊邏輯(fuzzylogic)小波(wavelet)、通信(communicati

9、on)、統計(statistics)等工具箱，這些工具箱給各個領域的研究和工程應用提供了有力的工具，借助于這些巨人肩上的工具”，各個層次的研究人員可直觀、方便地進行分析、計算及設計工作，從而大大地節(jié)省了時間。3.3.2. MATLAB研究數字水印的優(yōu)點集成了DCT、DWT等函數有豐富的小波函數和處理函數，這不僅方便了研究人員，而且使源程序簡潔明了、易實現。強大的數學運算功能。能夠方便、高效地實現音頻、視頻中的大量矩陣運算。提供了圖像處理工具箱、小波分析工具箱、數字信號處理工具箱。用來編制跨數字圖像處理技術、數字信號處理等多學科的數字水印技術是非常好的選擇。MATLAB與目前最強大的編程工具一一

10、VisualC+具有良好的接口。3.3.3. MATLAB函數介紹在介紹函數之前，我們必須明確一點：作水印程序時，處理的圖像數據是二維信號，而聲音信號是一維信號。這里，我們僅僅簡單介紹與水印有關的函數。數據輸入輸出函數imread()和imwrite():可以讀寫bmp,jpg/jpeg,tif/tiff,png,hdf,pcx,wxd格式文件。讀索引文件時，還可以得到相應的調色板數據。auread()、auwrite()、wavread()和wavwrite():可以方便地讀寫au和wav文件,并可控制其中的位及頻率。圖像顯示imshow()：顯示一幅圖像；imfinfo()：可以得到讀入圖

11、像的信息。如文件的大小、格式、格式版本號、圖像的高度、寬度、顏色類型(真彩色，灰度圖還是索引圖)等。變換頻函數對信號采用不同的變換，是實現頻域法水印的至關重要的一步，MATLAB中對一維信號和二維信號分別提供了各種變換和逆變換函數。離散余弦變換(DCT)DCT(離散余弦變換)；IDCT(逆向離散余弦變換)；DWT(離散小波變換)IDWT(離散小波變換)；dct(),dct2():分別實現一維信號和二維信號的idct(),idct2():分別實現一維信號和二維信號的離散小波變換(DWT)dwt(),dwt2():分別實現一維信號和二維信號的idwt(),idwt2():分別實現一維信號和二維信號

12、的Wavedec2():多級二維小波分解函數；Waveinfo():提供小波包中所有的小波信息；攻擊函數對算法進行攻擊測試是對水印魯棒性檢測的一種重要手段，一個好的水印算法必須經過各種攻擊測試才能對之做出客觀的評價。MATLAB中的許多函數可以直接用來做攻擊測試。旋轉：rotate()可以對圖像進行任意角度的旋轉；剪裁：imcrop()可以按精確定位的各點坐標進行剪裁；濾波：filter()和filter2()可實現對一維信號和二維信號的濾波；抖動：dither()對圖像進行抖動；抖動攻擊考驗水印魯棒性的一個很好的攻擊；jpeg壓縮：imwrite()中jpg和quality參數能對圖像進行可

13、控jpg壓縮；加各種噪聲：imnoise()可以對圖像加入各種噪聲，如白噪聲、椒鹽噪聲等，加入噪聲是對水印魯棒性考驗的一種常見的攻擊；放大/縮?。篿mresize()可以以指定的插值方法來對圖像進行放大和縮小。第四章：圖像數字水印技術的實現4.1基于離散余弦變法(DCT實現數字水印技術打開原始及水印圖像：subplot(2,2,1)I=uigetfile(*.bmp,打開原始彩色圖像文件,);RGB=imread(I);image(RGB);title(原始彩色圖像,);subplot(2,2,2)I=uigetfile(*.bmp,打開水印灰度圖像文件);imshow(I);title(灰度

14、水印圖像,);subplot(2,2,3)H=imread(I);J=dct2(H);imshow(log(abs(J),),colorbar;title(水印圖像經DCT變換后能量分布情況,)運行結果：水印全過程：0%水印加入程序Q=input(請輸入放縮因子的值(建議小于1):Q=)subplot(2,3,1)RGB=imread(浙江臺州學院,jpg);imshow(RGB);title(原始圖像);subplot(2,3,2)N=dct2(RGB(:,:,3);imshow(log(abs(N),),colorbar;title(Y分量能量分布,);subplot(2,3,4)I=im

15、read(lena1,bmp);imshow(I);title(灰度水印圖像,);subplot(2,3,5)M=dct2(I);imshow(log(abs(M),),colorbar;title(水印能量分布,);subplot(2,3,6)J=M(1:128,1:128);J(128:364,128:400)=0;J=rot90(J);J=rot90(J);J(365:600,401:750)=0;J=rot90(J);J=rot90(J);N=N+Q*J;K=idct2(N);RGB(:,:,3)=K;imshow(RGB);title(加入水印后圖像,);%水印提取程序subplot

16、(2,3,3)RGB1=imread(浙江臺州學院,jpg);N=dct2(RGB(:,:,3);M=dct2(RGB1(:,:,3);M=(N-M)/Q;B=idct2(M(236:365,350:401);Y=mat2gray(B);imshow(Y);title（提取的水印圖像,）運行結果：水印全過程(經剪切檢測水印)%水印加入程序Q=input(請輸入放縮因子的值(建議小于1):Q=)subplot(3,3,1)RGB=imread(MM,jpg);imshow(RGB);title(原始圖像);subplot(3,3,2)imshow(RGB(:,:,3);title(B分量);su

17、bplot(3,3,3)N=dct2(RGB(:,:,3);imshow(log(abs(N),),colorbar;title(B分量能量分布,);subplot(3,3,4)I=imread(lena1,bmp);imshow(I);title(灰度水印圖像,);subplot(3,3,5)M=dct2(I);imshow(log(abs(M),),colorbar;title(水印能量分布,);subplot(3,3,7)J=M(1:128,1:128);J(128:464,128:364)=0;J=rot90(J);J=rot90(J);J(465:800,365:600)=0;J=r

18、ot90(J);J=rot90(J);N=N+Q*J;K=idct2(N);RGB(:,:,3)=K;imshow(RGB);title(加入水印后圖像);subplot(3,3,8)I=imcrop(RGB,11598798);imshow(I);subplot(3,3,9)%水印提取程序subplot(3,3,6)RGB1=imread(MM,jpg);J=RGB1(:,:,3);X=J(1:799,1:599);N=dct2(I(:,:,3);M=dct2(X);M=(N-M)/Q;B=idct2(M(337:464,237:364);Y=mat2gray(B);imshow(Y);ti

19、tle(經放縮后提取的水印圖像)運行結果：水印全過程(經空域壓縮檢測水印)程序源代碼%水印加入程序Q=input(請輸入放縮因子的值(建議小于1):Q=)P=input(請輸入您所希望的圖像放縮系數值(建議取值不要小于0.5):P=)subplot(3,3,1)RGB=imread(浙江臺州學院,jpg);imshow(RGB);title(原始圖像);subplot(3,3,2)imshow(RGB(:,:,3);title(B分量);subplot(3,3,3)N=dct2(RGB(:,:,3);imshow(log(abs(N),),colorbar;title(B分量能量分布,);su

20、bplot(3,3,4)I=imread(lena1,bmp);imshow(I);title(灰度水印圖像,);subplot(3,3,5)M=dct2(I);imshow(log(abs(M),),colorbar;title(水印能量分布,);subplot(3,3,7)J=M(1:128,1:128);J(128:364,128:400)=0;J=rot90(J);J=rot90(J);J(365:600,401:750)=0;J=rot90(J);J=rot90(J);N=N+Q*J;K=idct2(N);RGB(:,:,3)=K;imshow(RGB);title(加入水印后圖像,

21、);subplot(3,3,8)I=imresize(RGB,Pnearest);imshow(I);title(壓縮P倍圖像)；subplot(3,3,9)J=imresize(I,1/P,nearest);imshow(J);title(再放大P倍還原圖像,)%水印提取程序subplot(3,3,6)RGB1=imread(浙江臺州學院,jpg);N=dct2(J(:,:,3);M=dct2(RGB1(:,:,3);M=(N-M)/Q;B=idct2(M(236:365,350:401);Y=mat2gray(B);imshow(Y);title(,經放縮后提取的水印圖像,):運行結果：4

22、.2圖像水印的dwt算法%以下是水印提取算法clearall;clc;%呆存時間start_time=cputime;figure(1);%讀出原始圖像subplot(1,2,1);input=imread(2.jpg);imshow(input);title(原始圖像)；%讀出水印圖像subplot(1,2,2);watermarked_image=imread(watermarked.bmp);imshow(watermarked_image,口)；title(水印圖像)；%三色分離input=double(input);watermarked_image=double(watermark

23、ed_image);inputr=input(:,:,1);watermarked_imager=watermarked_image(:,:,1);inputg=input(:,:,2);watermarked_imageg=watermarked_image(:,:,2);inputb=input(:,:,3);watermarked_imageb=watermarked_image(:,:,3);%水印圖像R勺分解Cwr,Swr=WAVEDEC2(watermarked_imager,2,haar);洞彳象R的分解Cr,Sr=WAVEDEC2(inputr,2,haar);%水印圖像G勺分

24、解Cwg,Swg=WAVEDEC2(watermarked_imageg,2,haar);洞彳象R的分解Cg,Sg=WAVEDEC2(inputg,2,haar);%水印圖像B的分解Cwb,Swb=WAVEDEC2(watermarked_imageb,2,haar);洞彳象B的分解Cb,Sb=WAVEDEC2(inputb,2,haar);%提取水印小波系數%提取水印R勺小波系數r=0.06;fork=0:3whr(k+1,:)=Cwr(1+size(Cwr,2)/4+k*size(Cwr,2)/16:size(Cwr,2)/4+(k+1)*size(Cwr,2)/16)-Cr(1+size

25、(Cr,2)/4+k*size(Cr,2)/16:.size(Cr,2)/4+(k+1)*size(Cr,2)/16);wvr(k+1,:)=Cwr(1+size(Cwr,2)/2+k*size(Cwr,2)/16:size(Cwr,2)/2+(k+1)*size(Cwr,2)/16)-Cr(1+size(Cr,2)/2+k*size(Cr,2)/16:.size(Cr,2)/2+(k+1)*size(Cr,2)/16);wdr(k+1,:)=Cwr(1+3*size(Cwr,2)/4+k*size(Cwr,2)/16:3*size(Cwr,2)/4+(k+1)*size(Cwr,2)/16)

26、-Cr(1+3*size(Cr,2)/4+k*size(Cr,2)/16:3*size(Cr,2)/4+(k+1)*size(Cr,2)/16);endwhr=(whr(1,:)+whr(2,:)+whr(3,:)+whr(4,:)/(4*r);wvr=(wvr(1,:)+wvr(2,:)+wvr(3,:)+wvr(4,:)/(4*r);wdr=(wdr(1,:)+wdr(2,:)+wdr(3,:)+wdr(4,:)/(4*r);war=(Cwr(1:size(Cwr,2)/16)-Cr(1:size(Cr,2)/16)/r;%提取水印小波系數%提取水印G勺小波系數g=0.03；fork=0:

27、3whg(k+1,:)=Cwg(1+size(Cwg,2)/4+k*size(Cwg,2)/16:size(Cwg,2)/4+(k+1)*size(Cwg,2)/16)-Cg(1+size(Cg,2)/4+k*size(Cg,2)/16:size(Cg,2)/4+(k+1)*size(Cg,2)/16);wvg(k+1,:)=Cwg(1+size(Cwg,2)/2+k*size(Cwg,2)/16:size(Cwg,2)/2+(k+1)*size(Cwg,2)/16)-Cg(1+size(Cg,2)/2+k*size(Cg,2)/16:.size(Cg,2)/2+(k+1)*size(Cg,2

28、)/16);wdg(k+1,:)=Cwg(1+3*size(Cwg,2)/4+k*size(Cwg,2)/16:3*size(Cwg,2)/4+(k+1)*size(Cwg,2)/16)-Cg(1+3*size(Cg,2)/4+k*size(Cg,2)/16:3*size(Cg,2)/4+(k+1)*size(Cg,2)/16);endwhg=(whg(1,:)+whg(2,:)+whg(3,:)+whg(4,:)/(4*g);wvg=(wvg(1,:)+wvg(2,:)+wvg(3,:)+wvg(4,:)/(4*g);wdg=(wdg(1,:)+wdg(2,:)+wdg(3,:)+wdg(4

29、)/(4*g);wag=(Cwg(1:size(Cwg,2)/16)-Cg(1:size(Cg,2)/16)/g;%提取水印小波系數%提取水印B的小波系數b=0.12;fork=0:3whb(k+1,:)=Cwb(1+size(Cwb,2)/4+k*size(Cwb,2)/16:size(Cwb,2)/4+(k+1)*size(Cwb,2)/16)-Cb(1+size(Cb,2)/4+k*size(Cb,2)/16:size(Cb,2)/4+(k+1)*size(Cb,2)/16);wvb(k+1,:)=Cwb(1+size(Cwb,2)/2+k*size(Cwb,2)/16:size(Cwb

30、,2)/2+(k+1)*size(Cwb,2)/16)-Cb(1+size(Cb,2)/2+k*size(Cb,2)/16:size(Cb,2)/2+(k+1)*size(Cb,2)/16);wdb(k+1,:)=Cwb(1+3*size(Cwb,2)/4+k*size(Cwb,2)/16:3*size(Cwb,2)/4+(k+1)*size(Cwb,2)/16)-Cb(1+3*size(Cb,2)/4+k*size(Cb,2)/16:3*size(Cb,2)/4+(k+1)*size(Cb,2)/16);endwhb=(whb(1,:)+whb(2,:)+whb(3,:)+whb(4,:)/