基於矩陣奇異值分解的水印算法
一.實驗目的
了解基於矩陣奇異值分解的圖像數字水印技術，掌握基於矩陣奇異值分解的圖像水印算法原理，設計並實現一種基於矩陣奇異值分解的數字水印算法。
二.實驗條件
(1) Windows 10或7操作系統；
(2) MATLAB 2014b；
(3)圖像文件
三.實驗原理
1.矩陣的奇異值分解（SVD）與圖像矩陣的能量
矩陣的奇異值分解變換是一種正交變換，它可以將矩陣對角化。我們知道任何一個矩陣都有它的奇異值分解，對於奇異值分解可用下面的定理來描述。
定理
設A是一個秩為 r的 
矩陣，則存在正交矩陣 和 ，使得 
，其中
，這里 
，
是矩陣
對應的正特征值。稱 
為A的奇異值分解， 
（ 
）稱為A 的奇異值。
矩陣的 F（Frobenius）范數定義為
， 由於
所以 
上式表明矩陣F 范數的平方等於矩陣的所有奇異值的平方和。對於一幅圖像，通常用圖像矩陣的 F范數來衡量圖像的能量，所以圖像的主要能量集中在矩陣哪些數值較大的奇異值上。
定理（Weyl定理）
設 A為一個大小為
的矩陣，
，
是矩陣A 的一個擾動，假設矩陣 A，B 的奇異值分別為 
和 ，
，
是矩陣 
的最大奇異值，則有
，其中 
表示2—范數。
由Weyl定理可知，當圖像被施加小的擾動時圖像矩陣的奇異值的變化不會超過擾動矩陣的最大奇異值，因此基於矩陣奇異值分解的數字水印算法具有很好的穩定性，能夠有效地抵御噪聲對水印信息的干擾。由上述矩陣奇異值分解的性質可知，圖像矩陣奇異值分解的穩定性非常好。當圖像加入小的擾動時，其矩陣奇異值的變換不超過擾動矩陣的最大奇異值。基於矩陣奇異值分解的數字水印算法正是將想要嵌入的水印信息嵌入到圖像矩陣的奇異值中，如果在嵌入水印的過程中選擇一個嵌入強度因子來控制水印信息嵌入的程度，那么當嵌入強度因子足夠小時，圖像在視覺上不會產生明顯的變化。
2 水印嵌入
設一副圖像對應的矩陣A 大小為
，需要嵌入的水印對應的矩陣 W大小為
，自然地有 
，
。在矩陣 A的左上角取一個大小為 
的子塊
。首先對 進行奇異值分解，得到
，其中 
是的奇異值矩陣。我們的目標就是將水印 
嵌入到矩陣
中，在這里定義一個描述水印嵌入過程的參數
，稱為嵌入強度因子，則水印嵌入的過程表示為
。可以看出，矩陣 
包含了所有的水印信息，水印信息的能量反映在 的奇異值當中。對進行奇異值分解，得到 
，則
反映了嵌入水印的圖像的全部信息，由此得到子塊 嵌入水印后的圖像子塊
，這樣就完成了水印的嵌入。
上述水印嵌入算法的基本過程可以表示為 
嵌入強度因子 的意義在上式中一目了然，它衡量了水印對原圖像的擾動情況。在水印嵌入時，選擇合適的嵌入強度因子是十分重要的。小的嵌入強度因子有利於水印的透明性，但嵌入的水印信息容易受到外界噪聲的干擾，如果噪聲強度足夠大，則可能使水印信息被噪聲淹沒而完全丟失，導致提取水印時無法得到水印的全部信息。大的嵌入強度因子有利於增強算法的魯棒性，即使在噪聲較強的情況下水印信息也不會受到很大的影響，但是過大的嵌入強度因子可能對原矩陣的奇異值產生較大的影響，有可能破壞水印的透明性，影響圖像的質量。因此，在水印嵌入時要選擇適當的嵌入強度因子使得水印圖像的不可覺察性與魯棒性達到最佳。
3 水印提取
水印提取是上述水印嵌入過程的逆過程。在水印提取時，假設我們得到的是受擾動的圖像矩陣
，首先對 
進行奇異值分解
由此得到包含有全部水印信息的奇異值矩陣 
，然后利用水印嵌入時的矩陣 
，得到
由上述水印嵌入的算法可得
這樣就得到了水印的信息
，其中 
為水印嵌入時所用的嵌入強度因子，
為原圖像 
的奇異值矩陣。
四.實驗步驟
1.打開matlab軟件，新建腳本文件並保存
2.輸入程序
3.編譯通過后點擊運行
五.實驗代碼
clc, clear A=imread('tu5.bmp'); %讀入載體文件 W=imread('tu6.bmp'); %讀入水印文件 [m1,m2,m3]=size(W); %給出矩陣W的維數 A0=A([1:m1],[1:m2],:); %在矩陣A的左上角選取與W同樣大小的子塊 A0=double(A0); W=double(W); %進行數據類型轉換 a=0.05; %嵌入強度因子為0.05 for i=1:3 [U1{i},S1{i},V1{i}]=svd(A0(:,:,i)); %對載體R，G，B層分別進行奇異值分解 A1(:,:,i)=S1{i}+a*W(:,:,i); %計算A1矩陣 [U2{i},S2{i},V2{i}]=svd(A1(:,:,i)); %對A1的各層進行奇異值分解 A2(:,:,i)=U1{i}*S2{i}*V1{i}'; %計算A2矩陣 end AW=A; %整體水印合成圖片初始化 AW([1:m1],[1:m2],:)=A2; %左上角替換成水印合成子塊，水印嵌入完成 AW=uint8(AW); W=uint8(W); %變換回原來的數據類型 subplot(1,3,1), imshow(A) %顯示載體圖片 subplot(1,3,2), imshow(W) %顯示水印圖片 subplot(1,3,3), imshow(AW) %顯示嵌入水印的合成圖片 %以下是水印的提出 AWstar=imnoise(AW,'gaussian',0,0.01); %加入高斯噪聲 A2star=AWstar([1:m1],[1:m2],:); %提出子塊 A2star=double(A2star); %進行數據類型轉換 for i=1:3 [U3{i},S2star{i},V3{i}]=svd(A2star(:,:,i)); %奇異值分解 A1star(:,:,i)=U2{i}*S2star{i}*V2{i}'; %計算A1* Wstar(:,:,i)=(A1star(:,:,i)-S1{i})/a; %計算W* end for i=1:3 Wstar(:,:,i)=medfilt2(Wstar(:,:,i)); %對提取水印的R,G,B層分別進行中值濾波 end Wstar=uint8(Wstar); %進行類型轉換 figure, subplot(1,2,1), imshow(AWstar) %顯示被噪聲污染的合成圖片 subplot(1,2,2), imshow(Wstar) %顯示提出的水印圖片
