图像数字水印+matlab程序
-
第三章
图像数字水印的方案
3.1
图像数字水印的技术方案
在数据库中存储在国际互联网上传输的水印图像一般会被压缩，
有时达到
很高的压
缩比。因此，数字水印算法所面临的第一个考验就是压缩。
JPEG
和
EZW
（
Embedded
Zero-Tree Wavelet
< br>）压缩是最常见的两种压缩方法。
JPEG
是基于离散余
弦变换域的压缩
方法，而
EZW
是基于
小波变换域的压缩方法。前人的研究证明采用与压缩算法相同的
变换域水印方法，
对于压缩的稳健性较强。
因此，
我研究图像文
件水印算法主要集中在
变换域算法及利用人眼视觉特性上。
数字水印的嵌入要求即要考虑视觉透明性，
< br>又要保证嵌入水印后图像的稳健性，
这
两个方面存在着矛
盾。
保证视觉透明性，
就要将水印嵌入到人眼不敏感区，
也就是嵌入
到图像的高频分量中。
而多数图像
处理方法对于图像高频部分的损坏程度较高，
如有损
压缩、
高频滤波等。
水印很容易在经历图像处理的过程中丢失。这样，
则无法保证图像
数字水印的稳健性。如果要获得很好的稳健性，数字
水印应加在人眼敏感的低频部分，
图像的大部分能量集中在低频部分，
< br>如果对于低频部分进行处理，
水印固然会失去，
而
图像也没有了利用价值，
然而，
水印的嵌入会
对图像的质量有非常大的影响，
这又无法
保证视觉透明性。
数字水印算法的实现基本分为三个部分
：
宿主图像的变换，
水印的嵌入和水印的检
测，分别描述如下。
3.2
基
于
DCT
域的图像数字水印技术
离散余弦变换
（
Discrete
Cosine Transform
）
属于正交变换图像编码方
法中的一种。
正交变换图像编码始于
1968
< br>年。当时安德鲁斯（
Andrews
）等人发现大多数自
然图像
的高频分量相对幅度较低，
可完全舍弃或者只用少数码字
编码，
提出不对图像本身编码，
只对其二维傅立叶
（
DFT
）
系数进行编码和
传输。
但
DFT
是一种正交变换，
运算量很大，
常常使实时处理发生困难，第二年他们就用
Walsh-Hadamard
变换（
WHT
）取代
DFT
可
以
使运算量明显减少，这是因为
WHT
变换只有加减法而无需乘法
。但是更有意义的是
离散余弦变换和离散正旋变换的出现，
它们
具有快速算法，
精确度高。
其中最重要的是
1974
年提出的
DCT
，因为其
变换矩阵的基向量很近似于托伯利兹矩阵的特征向量，而
托伯利兹矩阵又体现了人类语言
及图像信号的相关性。
因此，
DCT
常
常被认为是语音与
图像信号变换的准最佳变换。
图像是二维的，所以在研究时主要用到二维
DCT
，
以及二维
IDCT
来对图像进行处
理。
3.2.1
离散余弦变换（
DCT
）的定义
<
/p>
数字图像
X
（
m
,n
）是具有
M
行
N
列的一个矩阵。为了同时减弱或去除图像数据
相关性，可
以运用二维
DCT
，将图像从空间域转换到
DCT
变换域。
根据定义，二维
离散余弦变换（
DCT
）定义如下：
Y
(
k
,
l
)
2
MN
c
(
k
)
c
(
l
)
X
< br>(
m
,
n
)
cos
m
0
n
0
M
1
N
1
(
2
m
1
)
k
(
2
n
< br>
1
)
l
cos
2
M
2
N
式中：
m
，
k=0,1,…,M
-1; n
，
l
=0,1,…,N
-1
。
< br>
1
2
k
0
其中函数
c
(
k
)
k
1
,<
/p>
2
,
,
M
1
1
1
2
k
0
c
(
l
)
1
k
<
/p>
1
,
2
,
,
N
1
二维逆离散余弦变换（
IDCT
）的定义如下：
X
(
m
,
n
)
2
MN
M
1
< br>N
1
K
0
L
0
c
(
k<
/p>
)
c
(
l
)
Y
(
k
,
l
)
cos
(
2
m
1
)
k
(
2
n
1
)
l
<
/p>
cos
2
M
2<
/p>
N
式中：
m
，<
/p>
k=0,1,…,M
-1; n
，
l
=0,1,…N
-1
。
3.2.2
离散余弦变换的特点
在基于
DCT
的变换编码中，图像是先经分块（
8
×
8
或
16
×
16
）后再经
DC
T
，这种
变换是局部的，
只反映了图像
某一部分的信息。
当然也可以对整幅图像的特点，
但是运
算速度比分块
DCT
要慢。图像经
DCT
后，得到的
DCT
图像有三个特点：
一是系数值全部集中到
0
值附近（从直方图统计的意义上）
，动态范围很小，这说
明用较小的量化比特数即可表示
DCT
系数；
二是
DCT
< br>变换后图像能量集中在图像的低频部分，即
DCT
图像中
不为零的系数大
部分集中在一起（左上角）
，因此编码效率很高
。
三是没有保留原图像块的精细结构，从中反映不了原图像块
的边缘、轮廓等信息，
这一特点是由
DCT
缺乏时局域性造成的。
如下左图
3
—
1
是原始图像经过
DCT
变换后的系数图像为图
3
—
2
。
两条线划分出图
像的低频、中频和高频分别所在的矩形区域。可以看出，图像
DCT
变换后大部分参数
接近于零，
只有左上角的低频
部分有较大的数值，
中频部分参数值相对较小，
而大部分
高频参数值非常小，接近于零。
<
/p>
图
3
—
1
：原图像
图
3
—
2
：变换后的系数图像
3.2.3
离散余弦变换的数字水印算法
根据离
散余弦变换后的参数性质，本文采用了以
ZigZag
方式重排
变换域系数的方
法，选出中频分量，
用数字水印序列对其进行非
线性调制。水印检测时，待检测图像仍
按比方式选择变换域系数，与待水印进行相关运算
，与阈值比较来判断是否所含水印。
离散余弦域的数字水印算
法的具体实现分为三步：
宿主图像的变换，
数字水印的嵌
入，数字水印的检测。
3.2.3.1 <
/p>
宿主图像的
DCT
变换
< br>
对于
N
×
N
大小的
256
灰度级的宿主图像
I
进行
N
×<
/p>
N
二维离散余弦变换（
DCT
）
。
以
ZigZag
方式对于
DCT
变换后的图像频率系数重新
排列成一维向量
Y={y
1
, y
2
,…y
N×N
}
.
并取出序列中第
L+1
到
L+M
的中频系数部分
,
得到
Y
L
={
Y
L+1
,
Y
L+2
,…,
Y
L+M
}
3.2.3.2
数字水印的嵌入
<
/p>
假设数字水印
W
为一服从标准正态分布的
随机实数序列，用数字序列表示为
W={W
1
< br>，
W
2
，
…W
M
}
。用
W
对
Y
序
列中第
L+1
到
L+M
的中频系数部分的值进行修改，
按以下公式进行：
经过修改的系数序列
Y
′
={
Y
1
′
, <
/p>
Y
2
′
,...
Y
′
N×N
}
以
ZigZag
逆变换形式重组，再<
/p>
进行
N
×
N D
CT
逆变换，得到嵌有数字水印的图像
I
′。
y
i
'
y
i<
/p>
'
2
y
i
y
i
y
i
w
i
i
L
或
i
L
m
L<
/p>
i
L
M
,
0
3.2.3.3
数字水印的检测
待检测的可能含有水
印的图像
I
。假设
I
未损失大量信息，可以近似认为
I
= I
′。
在此假设下可以运用统计的方法来检测水印。
（
1
）待检水印域待检图像中频系数相关性的测定
<
/p>
同样对
I
′进行
DCT
变换，
以
ZigZag
方式将
DCT
系数排成一维向量
Y
{
Y
< br>1
Y
2
Y
N×N
。由
于假设
I
=I
′，则
Y
=
Y
′。
取出
Y
（
等于
Y
′）
中第
L+1
到
L+M
的中频系数部分
Y
L
Y
L+1
Y
L+2
,...
Y
L+M
。
假设待检测的数字水印
X={X
1
，
X
2
< br>，
... X
M
}
为一符合标准正态分布的实数伪随机序列。
则可以通过待检水印与图像中频系
数作相关运算来判断是否所加入了水印。
只有在待检
水印为所加
入的水印时，才能得到较大的相关值。否则相关值很小，接近于零。
< br>’
1
Z
cov(
Y
L
'
,
X
)
M
(
Y
L<
/p>
i
X
i
)
i
1
M
1
M
((
Y
L
< br>
i
X
i
)
(
Y
L
2
i
W
i
Xi
))
用符号
E
表示数
学期望，得到：
（
2
）阈值的确定
E
(
z
)
E
(
< br>Y
2
)
X
W
X
W
E
(
z
)
0
E
(
z
)
0
没有水印存在
根据中心极限定理，
参照水印匹配与不匹配两种情况得到阈值为
T
z
E
(
Y
L
2
)
/
2
。
由于原始图像难以得到，因
此从实用性出发，阈值定义为：
T
Z
E
(
Y
)
/
2
2
L
M
Y
i<
/p>
1
M
L
i
综上所述，满
足
Z
T
Z<
/p>
或
匹配的水印。
2
0
.
5
时，则表明检测到匹配水印。否则，未检测到
2
T
Z
3.3
MATLAB
工具简介
3.3.1.
简介
Matlab
< br>是当前在国内外十分流行的工程设计和系统仿真软件包。它是
MathWorks
公司于
1982
年推出的一套高性能的
数值计算和可视化软件，
它集数值分析、
矩阵运算、
信号处理和图形显示于一体，构成了一人方便的、界面友好的用户环境。
Matlab
的推出得到了各个领域专家、学者的广泛关注，其
强大的扩展功能为各个领域
的应用提供了基础。
由各个专家学者
相继推出了
MATLAB
工具箱，
其中
的信号处理
(signal
processing)
、控制系统
(control
system)
、神经网络
(neural
netw
ork)
、图像处理
(image
p
rocessing)
、鲁棒控制
(robust
control)
、非线性系统控制设计
(non
linear
system
control
design)
、系统辨识
(system ide
ntification)
、最优化
(optimizatio
n)
、模糊逻辑
(fuzzy logic)
< br>、
小波
(wavelet)
、通
信
(communication)
、统计
(statistics)
等工具箱，这些工具箱给各个领域
的研究和工程应用提供了有力的工具，
借助于这些
“
巨人肩上的工具
”
，
各
个层次的研究
人员可直观、方便地进行分析、计算及设计工作，从而大大地节省了时间。
3.3.2.
MATLAB
研究数字水印的优点
①
集成了
D
CT
、
DWT
等函数有丰富的小波函数
和处理函数，这不仅方便了研究人员，
而且使源程序简洁明了、易实现。
②
强大的数学运算功能。
能够方便、高效地实现音频、视频中的大量矩阵运算。
③
提供了图像处理工具箱、小波分析
工具箱、数字信号处理工具箱。用来编制跨数字
图像处理技术、数字信号处理等多学科的
数字水印技术是非常好的选择。
④
MATLAB
与目前最强大的编程工具
——
Visual
C++
具有良好的接口。
3.3.3. MATLAB
函数介绍
在介绍函数之前，
我们必须明确一点
：
作水印程序时，
处理的图像数据是二维信号，
①
数据输入输出函数
imread
（）和
imwrite
（）
：可以读写
bmp,jpg/jpeg,
tif/tiff, png, hdf, pcx, wxd
格式文件。
读索引文件时，还可以得到相应的调色板数据。
a
uread
（）
、
auwrite
（）
、
wavread
（）
和
wavwrite
（）
：
可以方便地读写
au
和
wav
文件，
并可控制其中的位
及频率。
②
图像显示
imshow()
：
显示一幅图像；
i
mfinfo()
：
可以得到读入图像的信息。
如文件的大小、
格式、
格式版本号、图像的高度、宽度
、颜色类型（真彩色，灰度图还是索引图）等。
③
变换频函数
对信号采用不同的变换，
是实现频域法水印的至关重要的一步，
MATLAB
中对一维<
/p>
信号和二维信号分别提供了各种变换和逆变换函数。
离散余弦变换（
< br>DCT
）
dct(),dct
2():
分别实现一维信号和二维信号的
DCT
（离散余弦变换）
；
idc
t(),idct2():
分别实现一维信号和二维信号的
ID
CT
（逆向离散余弦变换）
；
离散小波变换（
< br>DWT
）
dwt(),dwt
2():
分别实现一维信号和二维信号的
DWT
（离散小波变换）
idwt(),idwt2():
分别实现一维信号和二维信号的
IDWT
（离散小波变换）
；
Wavede
c2():
多级二维小波分解函数；
Waveinfo():
提供小波包中所有的小波信息；
④攻击函数
对算法进行攻击测试是对
水印鲁棒性检测的一种重要手段，
一个好的水印算法必须经过
各
种攻击测试才能对之做出客观的评价。
MATLAB
中的许多函
数可以直接用来做攻击测
试。
旋转：
rotate
（）可以对图像进行任意角度的旋转；
剪裁：
imcrop()
可以按精确定位的各点坐标进行剪裁；
滤波：
filter()
和
filter2()
可实现对一维信号和二维信号的滤波；
抖动：<
/p>
dither()
对图像进行抖动；抖动攻击考验水印鲁棒性的一
个很好的攻击；
jpeg
压缩：
imwrite()
中
jpg
和
quality
参数能对图像进行可控
jpg
压缩；
加各种噪
声：
imnoise()
可以对图像加入各种噪声，如白噪声、
椒盐噪声等，加入噪声是
对水印鲁棒性考验的一种常见的攻击；
放大
/
缩小：
imresize()
可以以指定的插值方法来对图像进行放大和缩小。
而声音信号是一维信号。这里，我们仅仅简单介绍与水印有关的函数。
第四章：
图像数字水印技术的实现
4.1 <
/p>
基于离散余弦变法（
DCT
）实现数字水
印技术
①
打开原始及水印图像
：
subplot(2,2,1)
I=uigetfile('*.bmp','
打开原始彩色图像文件
');
RGB=imread(I);
image(RGB);
title('
原始彩色图像
');
subplot(2,2,2)
I=uigetfile('
*.bmp','
打开水印灰度图像文件
');
imshow(I);
title('
灰度水印图像
');
subplot(2,2,3)
H=imread(I);
J=dct2(H);
imshow(log(abs(J)),[]),colorbar;
title('
水印图像经
DCT
< br>变换后能量分布情况
')
运行结果：
②
水印全过程：
0%
水印加入程序
< br>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=imread('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(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
分量
');
subplot(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)