CN102156955A - 基于直方图邻域的鲁棒可逆水印嵌入与提取方法 - Google Patents

Abstract

本发明公开了一种基于直方图邻域的鲁棒可逆水印嵌入与提取方法，主要解决现有方法综合性能差的缺陷。该方法首先将原始图像分块，利用嵌入级别筛选出一组候选块；然后计算候选块的特征值，确定峰值对，并通过峰值对邻域选择的方式获得一组可靠的嵌入块；最后，通过修改嵌入块的特征值实现水印嵌入；提取时，获得待检测图像的嵌入块，利用聚类算法将嵌入块的特征直方图值序列进行分类，得到一组分类序列，然后，通过比较嵌入块的特征值与分类序列来提取水印，随后，修改待检测图像中嵌入块的特征值获得重构图像。本发明具有对JPEG和JPEG2000两种有损压缩攻击可靠的鲁棒性，同时容量大、不可见性好，整体性能优于同类方法，可用于数字图像的版权保护。

Description

基于直方图邻域的鲁棒可逆水印嵌入与提取方法

技术领域

本发明属于多媒体信息安全领域，涉及数字图像水印嵌入与提取方法，可用于网络环境中数字图像的内容认证、版权保护、侵权追溯。

背景技术

随着计算机及网络通信技术的飞速发展，以多媒体信息为代表的数字技术已经广泛渗透到军事、工业、医疗、教育等各个领域，这标志着人类社会进入了一个全新的数字化信息时代。人们借助于扫描仪、打印机等先进的电子设备和成熟的计算机软件可以方便、迅捷地将图像、视频等数字信息传送到世界的各个角落。一方面，不仅提高了信息表达的效率、准确性和便捷性，而且实现了信息存取、发布和传播的“数字化”与“网络化”。然而，另一方面，这种便利又会给有恶意的个人或团体带来可乘之机，使他们在没有得到作品所有者授权的情况下任意复制、修改、传播有版权的内容，甚至非法用作商业或军事用途，这不仅给版权所有者造成巨大的经济损失，而且对国家安全带来严峻的挑战。因此，如何在网络环境下对数字媒体实施版权保护、内容认证、侵权追溯等已经成为迫在眉睫的信息安全问题。

针对上面的问题，数字水印技术作为一项富有潜力的解决手段应运而生，并成为近年来工业界和学术界共同关注的热点。数字水印技术是利用数字内嵌的方法把一些秘密信息隐藏在图像、声音、文档、视频等数字媒体中。通常称嵌入的秘密信息为数字水印，它可以是版权标志、用户序列号或者产品认证信息等。通过对数字水印的检测与分析，以识别多媒体产品的所有者或安全标识信息，起到版权保护、内容认证和侵权追溯的作用。然而，大多数情况下，水印数据在嵌入到数字媒体的过程中都会给原始媒体造成一些不可逆的失真。即使这些失真微乎其微，不易察觉，但在医学影像、法庭证据照片、电子票据、军事及遥感图像等领域几乎是不允许的。因此，在满足版权保护与内容认证的同时，希望在水印数据提取之后能保持原始媒体内容的完整性。此外，由于多媒体产品在传输过程中往往会存在信道噪声和干扰，如有损压缩的影响，嵌入的水印又需要具有抵抗攻击的鲁棒性。这种水印方案称为鲁棒可逆水印，它以其特有的优势引起了广大研究者们的浓厚兴趣，成为国际信息安全领域新兴的研究领域。

根据水印嵌入模型的不同，现有的鲁棒可逆图像水印方法可以分为两类，见文献“AnL L，Gao X B，Deng C，and Ji F.Robust lossless data hiding：analysis and evaluation.In Proc.International Conference on High Performance Computing & Simulation，Caen，France，June 28-July 2，2010”。

第一类是基于直方图旋转的方法，该方法首先将原始图像分成大小相同且互不重叠的图像块，每个块中的像素又被随机地划分为两个相同大小的子区域并将其灰度直方图映射到圆上；最后通过旋转圆心到子区域质心的向量来嵌入水印。在接收端，根据两个子区域向量的差异来提取水印，并逆向操作恢复原始图像，见文献“De Vleeschouwer C，Delaigle J，and Macq B.Circular interpretation of bijective transformations in lossless watermarking formedia asset management.IEEE Trans.Multimedia，5(1)：97-105，2003”。虽然该方法对联合图像专家组JPEG压缩具有鲁棒性，但由于采用模256运算来避免像素溢出，导致含水印图像中存在大量的“椒盐”噪声，大大降低了含水印图像的视觉质量，并影响了该方法的鲁棒性。

第二类是基于直方图分布约束的方法，见文献“Zou D K，Shi Y Q，Ni Z C，and Su W.Asemi-fragile lossless digital watermarking scheme based on integer wavelet transform.IEEETrans.Circuits and Systems for Video Technology，16(10)：1294-1300，2006”和“Ni Z C，Shi YQ，Ansari N，Su W，Sun Q B，and Lin X.Robust lossless image data hiding designed forsemi-fragile image authentication.IEEE Trans.Circuits and Systems for Video Technology，18(4)：497-509，2008”。此类方法根据图像的直方图分布信息和待嵌入的水印有约束地修改图像的统计特性来实现水印嵌入。与第一类方法相比，该类方法虽然消除了“椒盐”噪声，提高了含水印图像的视觉质量，但是，由于误差纠错码的使用，使其仍存在容量低、可逆性与鲁棒性不可靠的问题。虽然文献“Gao X B，An L L，Li X L，and Tao D C.Reversibilityimproved lossless data hiding.Signal Processing，89(10)：2053-2065，2009”改进了Ni等人所提方法的可逆性，但还是没有最终彻底解决该类方法容量低和鲁棒性不可靠的问题。

发明内容

本发明的目的在于针对上述已有技术存在的不足，提出一种基于直方图邻域的鲁棒可逆水印嵌入与提取方法，以改善含水印图像的视觉质量，提高水印嵌入容量，增强鲁棒性的可靠程度，进而提高鲁棒可逆水印嵌入与提取方法的综合性能。

为实现上述目的，本发明的技术方案包括水印嵌入和水印提取。

一、基于直方图邻域的鲁棒可逆水印嵌入方法，包括如下步骤：

(1)将原始图像I分成大小相同且互不重叠的块，丢弃不能被完全划分的块，得到一组原始图像块序列B；

(2)将原始图像块序列B进行筛选，得到候选块序列A；

(3)计算候选块序列A的特征值并生成特征直方图，在特征直方图中确定峰值对W；

(4)用峰值对W在候选块序列A中查找嵌入块，得到一组嵌入块序列E以及嵌入块序列E在原始图像块序列B中的位置信息S；

(5)确定嵌入级别数L，L为正整数，并利用伪随机数产生函数生成二值的伪随机水印序列P；

(6)用嵌入级别数L修改嵌入块序列E的特征值，将水印序列P嵌入到嵌入块序列E中，得到含水印图像I^E：

(6a)从水印序列P中选择第j位水印p_j，p_j∈{0，1}，1≤j≤f，f表示水印序列P的长度；

(6b)修改嵌入块序列E中第i个块的特征值ξ_i，将水印p_j进行嵌入，嵌入规则为：

式中，

是嵌入水印后嵌入块序列E中第i个块的特征值，ξ_i是嵌入水印前嵌入块序列E中第i个块的特征值，1≤i≤σ，σ为嵌入块序列E中块的个数，p_j是水印序列中的第j位水印，L是嵌入级别数，ξ_l是峰值对W的左值，ξ_r是峰值对W的右值；

(6c)重复步骤(6a)～(6b)直到将水印序列P中的所有位水印都嵌入到原始图像I中，得到含水印图像I^E。

二、基于直方图邻域的鲁棒可逆水印提取方法，包括如下步骤：

1)将待检测图像I^S分成大小相同且互不重叠的块，得到一组待检测图像块序列C；

2)对待检测图像块序列C进行筛选，得到一组待检测嵌入块序列E^S；

3)计算待检测嵌入块序列E^S的特征值并生成特征直方图，将特征直方图中的值序列进行排序，得到一组有序的特征直方图值序列H；

4)利用K-means聚类算法将有序的特征直方图值序列H进行划分，得到一组分类序列Ω：

4a)利用水印嵌入过程中得到的峰值对W的左频率F_l和右频率F_r的取值，对K-means聚类算法进行初始化，初始化的值包括聚类个数K、聚类集合Ω₁，L，Ω_K，聚类中心D₁，L，D_K和迭代次数θ；

4b)按照特征直方图值序列H中每个值与K个聚类中心D₁，L，D_K之间的距离，将每个值分配到最相近的一个聚类集合中；

4c)更新K个聚类中心和迭代次数θ；

4d)重复步骤4b)～4c)直到找到满足条件

的聚类集合Ω₁，L，Ω_K或者使迭代次数θ达到100次，其中，h_i是特征直方图值序列H中的第i个值，1≤i≤μ，μ表示特征序列

中取不同值的元素个数，Ω_j是第j个聚类集合，

是第θ+1次迭代时第j个聚类集合的聚类中心，1≤j≤K，K表示聚类个数，将聚类集合Ω₁，L，Ω_K组成一个集合，得到分类序列Ω＝{Ω₁，L，Ω_K}，Ω_K是分类序列Ω中的第K个值；

5)将待检测嵌入块序列E^S的特征值与分类序列Ω进行比较，提取水印序列R，并用嵌入级别数L修改待检测嵌入块序列E^S的特征值，获得重构图像I^R。

本发明具有以下优点：

(1)本发明由于将原始图像分成大小相同且互不重叠的块，并在每个块中嵌入水印，使得水印容量大且易于控制；

(2)本发明由于将原始图像块序列进行筛选，避免了水印嵌入过程中的像素溢出，去除了含水印图像中的“椒盐”噪声，改善了含水印图像的视觉质量；

(3)本发明由于采用K-means聚类算法将有序的特征直方图值序列进行划分，并用待检测嵌入块序列的特征值与分类序列进行比较来提取水印，增强了鲁棒性的可靠程度；

(4)本发明由于采用阈值与峰值对来构造嵌入块序列进行嵌入，便于有效折衷水印容量与含水印图像的视觉质量；

(5)本发明由于采用嵌入级别数来控制水印嵌入强度，便于有效折衷鲁棒性与含水印图像的视觉质量；

(6)本发明由于改善了含水印图像的视觉质量，提高了水印嵌入容量，增强了鲁棒性的可靠程度，进而提高了鲁棒可逆水印嵌入与提取方法的综合性能。

附图说明

图1是本发明的水印嵌入流程框图；

图2是本发明的水印提取流程框图；

图3是本发明与现有一种方法生成的含水印图像的对比结果图。

具体实施方式

参照图1，本发明基于直方图邻域的鲁棒可逆水印嵌入方法，其实现步骤如下：

步骤1：将原始图像I分成大小相同且互不重叠的块，丢弃不能被完全划分的块，得到一组原始图像块序列B。

将大小为M×N的原始图像I按照m×n的块大小分成大小相同且互不重叠的块，其中n为偶数，当原始图像的长度M不能被分块的长度m整除或原始图像的宽度N不能被分块的宽度n整除时，在原始图像的边界处会产生块大小小于m×n的不能被完全划分的块，丢弃这些不能被完全划分的块，得到一组原始图像块序列B。

步骤2：将原始图像块序列B进行筛选，得到候选块序列A。

2.1)将原始图像块序列B中第i个块的所有像素分成两个大小相同的像素集合

和

1≤i≤v，v为原始图像块序列B中块的个数：

$Z_a^i=\left\{Z^{(x,y)}\right|Z^{(x,y)}\∈Z,\mathrm{mod}(x+y,2)=0\}$

$Z_b^i=\left\{Z^{(x,y)}\right|Z^{(x,y)}\∈Z,\mathrm{mod}(x+y,2)\≠0\},$

式中，Z^(x，y)是第i个块中位于第x行、第y列的像素的灰度值，1≤x≤m，1≤y≤n，m×n是块大小，Z是第i个块中所有像素灰度值的集合，mod(·)是取余数运算，mod(x+y，2)是计算x+y除以2所得到的余数；

2.2)从得到的两个像素集合

和

中分别按顺序选取一个像素，用这两个像素的灰度值构成一个像素对

其中

代表从

中所选择的第k个像素的灰度值，

代表从

中所选择的第k个像素的灰度值，1≤k≤m×n/2；

2.3)按如下方式对第i个块中的所有像素对

1≤k≤m×n/2进行判定：

$\left\{\begin{array}{c}{Z}_a^{(i,k)}<L,\\ Z_b^{(i,k)}>2^t-1-L\end{array}\right.$ 或者 $\left\{\begin{array}{c}{Z}_a^{(i,k)}>2^t-1-L,\\ Z_b^{(i,k)}<L\end{array}\right.,$

式中，L是嵌入级别数，t为表示原始图像时所用的二进制位数，取8；若第i个块中的所有像素对

1≤k≤m×n/2都不满足上述条件，认为第i个块为候选块；

2.4)重复步骤2.1)～2.3)直到将原始图像块序列B中每个块都判断完为止，得到一组候选块序列A。

步骤3：计算候选块序列A的特征值并生成特征直方图，在特征直方图中确定峰值对W候选块序列A的特征值并生成特征直方图，在特征直方图中确定峰值对W。

3.1)计算候选块序列A中每个块的特征值，组成特征序列Ψ＝(ξ₁，Kξ_i，Kξ_ω)，其中，ξ_i是第i个块的特征值，1≤i≤ω，ω为候选块序列A中块的个数，ξ_i表示为：

${\mathrm{\&xi;}}_i=\frac{2}{m\&times;n}\underset{k=1}{\overset{m\&times;n/2}{\mathrm{\&Sigma;}}}(Z_a^{(i,k)}-Z_b^{(i,k)})$

式中，m×n为块大小，

为从

中所选择的第k个像素的灰度值，

为从

中所选择的第k个像素的灰度值，1≤k≤m×n/2；

3.2)将特征序列Ψ＝(ξ₁，Kξ_i，Kξ_ω)中取值相同的元素合并，得到特征序列Ψ的值序列X，表示为X＝{x₁，L x_i，L x_η}，其中，x_i是值序列X中的第i个值，1≤i≤η，η表示特征序列Ψ中取不同值的元素个数；

3.3)计算值序列X中每个值在特征值序列Ψ中的频率：

F_i＝count({ξ_j|ξ_j＝x_i})

式中，F_i是值序列X中的第i个值x_i在特征值序列Ψ中的频率，该频率等于特征序列Ψ中取值为x_i的元素个数，count({ξ_j|ξ_j＝x_i})是计算集合{ξ_j|ξ_j＝x_i}中元素的个数，ξ_j是特征序列Ψ中的第j个元素，1≤j≤ω，ω为候选块序列A中块的个数，x_i是值序列X中的第i个值，1≤i≤η，η表示特征序列Ψ中取不同值的元素个数；

通过步骤3.2)和3.3)就可以生成候选块序列A的特征直方图，该直方图包括值序列X＝{x₁，L x_i，L x_η}和频率序列F＝{F₁，L F_i，L F_η}，其中，x_i是值序列X中的第i个值，F_i是频率序列F中的第i个值，1≤i≤η，η表示特征序列Ψ中取不同值的元素个数；

3.4)在特征直方图的频率序列F中找出当频率取最大值时所对应的值序列X中的特征值x_r，表示为：

$x_r=\underset{x_i}{\arg \max }{F}_i,$

式中，x_i是值序列X中的第i个值，F_i是频率序列F中的第i个值，1≤i≤η，η表示特征序列Ψ中取不同值的元素个数；

3.5)在特征直方图的频率序列F中找出当频率取次最大值时所对应的值序列X中的特征值x_l，表示为：

$x_l=\underset{x_i}{\arg \max }{F}_i,$ x_i≠x_r，

式中，x_i是值序列X中的第i个值，F_i是频率序列F中的第i个值，1≤i≤η，η表示特征序列Ψ中取不同值的元素个数，x_r是频率序列F中当频率取最大值时所对应的值序列X中的特征值；

3.6)根据特征值x_r和x_l生成峰值对W，表示为W＝[(ξ_l，F_l)，(ξ_r，F_r)]，其中ξ_l＝min(x_r，x_l)是峰值对W的左值，ξ_r＝max(x_r，x_l)是峰值对W的右值，x_r是频率序列F中当频率取最大值时所对应的值序列X中的特征值，x_l是频率序列F中当频率取次最大值时所对应的值序列X中的特征值，F_l是特征直方图的频率序列F中的最大值，表示峰值对W的左频率，F_r是特征直方图的频率序列F中的次最大值，表示峰值对W的右频率。

步骤4：用峰值对W在候选块序列A中查找嵌入块，得到一组嵌入块序列E以及嵌入块序列E在原始图像块序列B中的位置信息S。

4.1)生成一个向量S＝(s₁，K s_i，K s_v)，其中，s_i是向量S中的第i个值，1≤i≤v，v是原始图像块序列B中块的个数，将向量中的每个值初始化为0；

4.2)计算候选块序列A中第i个块的特征值ξ_i与峰值对W的左值ξ_l之间的距离d(ξ_i，ξ_l)＝|ξ_i-ξ_l|，以及特征值ξ_i与峰值对W的右值ξ_r之间的距离d(ξ_i，ξ_r)＝|ξ_i-ξ_r|，1≤i≤ω，ω为候选块序列A中块的个数；

4.3)当ξ_i≤ξ_l且d(ξ_i，ξ_l)小于或等于给定的阈值T或者当ξ_i≥ξ_r且d(ξ_i，ξ_r)小于或等于给定的阈值T时，认为第i块是嵌入块，记录该块在原始图像块序列B中的位置，即将向量S＝(s₁，K s_i，K s_v)中的第i个值s_i改为1；

4.4)重复步骤4.2)～4.3)，直到候选块序列A中所有块都判断完为止，得到一组嵌入块序列E以及嵌入块序列E在原始图像块序列B中的位置信息S。

步骤5：确定嵌入级别数L，L为正整数，并利用伪随机数产生函数生成二值的伪随机水印序列P。

首先，选择一个正整数作为嵌入级别数L；然后利用伪随机数产生函数rand(1，f)生成一个长度为f的二值伪随机数序列P＝{p₁，Kp_j，Kp_f}，将该序列作为水印序列，其中，p_j∈{0，1}，1≤j≤f，f表示水印序列P的长度，即水印序列中水印的位数，f小于或者等于嵌入块序列E中块的个数；

步骤6：用嵌入级别数L修改嵌入块序列E的特征值，将水印序列P嵌入到嵌入块序列E中，得到含水印图像I^E。

6.1}从水印序列P中选择第j位水印p_j，p_j∈{0，1}，1≤j≤f，f表示水印序列P的长度；

6.2)修改嵌入块序列E中第i个块的特征值ξ_i，将水印p_j进行嵌入，嵌入规则为：

式中，是嵌入水印后嵌入块序列E中第i个块的特征值，ξ_i是嵌入水印前嵌入块序列E中第i个块的特征值，1≤i≤σ，σ为嵌入块序列E中块的个数，p_j是水印序列中的第j位水印，L是嵌入级别数，ξ_l是峰值对W的左值，ξ_r是峰值对W的右值；其中，根据水印p_j的不同，修改ξ_i的具体实现如下：

6.2a)当p_j＝0时，将第i个块中所有像素的灰度值保持不变，使

6.2b)当p_j＝1时，如果ξ_i≥ξ_r，修改规则为：

式中，

是嵌入水印后第i个块的像素集合

中第k个像素的灰度值，是嵌入水印前第i个块的像素集合

中第k个像素的灰度值，L是嵌入级别数，

是嵌入水印后第i个块的像素集合

中第k个像素的灰度值，

是嵌入水印前第i个块的像素集合

中第k个像素的灰度值，t为表示原始图像时所用的二进制位数，取8；

6.2c)当p_j＝1时，如果ξ_i≤ξ_l，修改规则为：

6.3)重复步骤6.1)～6.2)直到将水印序列P中的所有位水印都嵌入到嵌入块序列E中，并根据嵌入块序列E在原始图像块序列B中的位置信息S，将嵌入水印的块和未嵌入水印的块进行组合，得到含水印图像I^E。

通过上述步骤1～步骤6可以实现水印的嵌入，将二值的伪随机水印序列P嵌入到原始图像中I，得到含水印图像I^E。

参照图2，本发明基于直方图邻域的鲁棒可逆水印提取方法，其实现步骤如下：

步骤A：将待检测图像I^S分成大小相同且互不重叠的块，得到一组待检测图像块序列C。

步骤B：对待检测图像块序列C进行筛选，得到一组待检测嵌入块序列E^S。

A1)利用水印嵌入过程中得到的嵌入块序列在原始图像块序列中的位置信息S＝(s₁，K s_i，K s_v)，s_i是向量S中的第i个值，s_i∈{0，1}，对待检测图像块序列C中的第i个块进行判定：如果s_i＝1，认为该块是待检测嵌入块，如果s_i＝0，认为该块不是待检测嵌入块，其中，1≤i≤v，v为原始图像块序列中块的个数，该个数与待检测图像块序列C中的块个数相同；

A2)重复步骤A1)直到待检测图像块序列C中的所有块都判断完为止，得到一组待检测嵌入块序列E^S。

步骤C：计算待检测嵌入块序列E^S的特征值并生成特征直方图，将特征直方图中的值序列进行排序，得到一组有序的特征直方图值序列H。

C1)计算待检测嵌入块序列E^S中每个块的特征值，组成特征序列

其中，

是第i个块的特征值，1≤i≤σ′，σ′为待检测嵌入块序列E^S中块的个数，

表示为：

式中，m×n为块大小，

为从第i个块的像素集合

中所选择的第k个像素的灰度值，

为从第i个块的像素集合中所选择的第k个像素的灰度值，1≤k≤m×n/2；

C2)将特征序列中取值相同的元素合并，得到特征序列

的值序列

表示为其中，

是值序列

中第i个值，1≤i≤μ，μ表示特征序列

中取不同值的元素个数；

C3)计算值序列中每个值在特征序列

中的频率：

式中，

是值序列

中的第i个值

在特征序列

中的频率，该频率等于特征序列

中取值为的元素个数，

是计算集合

中元素的个数，

是特征序列

中的第j个元素，1≤j≤σ′，σ′为待检测嵌入块序列E^S中块的个数，

是值序列

中的第i个值，1≤i≤μ，μ表示特征序列

中取不同值的元素个数；

通过步骤C2)和C3)就可以生成待检测嵌入块序列E^S的特征直方图，该直方图包括值序列，和频率序列

其中，

是值序列

中的第i个值，

是频率序列

中的第i个值，1≤i≤μ，μ表示特征序列

中取不同值的元素个数；

C4)将值序列

按照取值进行递增排序，得到一组有序的特征直方图值序列H＝[h₁，L h_i，L h_μ]，其中，h_i是有序的特征直方图值序列H中的第i个值，1≤i≤μ，μ表示特征序列中取不同值的元素个数，h₁＜L＜h_i＜L＜h_μ。

步骤D：利用K-means聚类算法将有序的特征直方图值序列H进行划分，得到一组分类序列Ω。

D1)利用水印嵌入过程中得到的峰值对W的左频率F_l和右频率F_r的取值，对K-means聚类算法的聚类个数K进行初始化：当F_l≠0且F_r≠0时，K取3；否则，K取2；

D2)对K-means聚类算法的聚类中心D₁，L，D_K进行初始化，当K＝3时，3个聚类中心分别为D₂＝0，

当K＝2时，2个聚类中心分别为

其中，h₁是有序的特征直方图值序列H中的最小值，

是向下取整运算，h_μ是有序的特征直方图值序列H中的最大值；

D3)将K个聚类集合Ω₁，L，Ω_K分别初始化为空集，使迭代次数θ取1；

D4)按照有序的特征直方图值序列H中每个值与K个聚类中心D₁，L，D_K之间的距离，将每个值分配到最相近的一个聚类集合中：

h_i∈Ω_j，如果

l＝1，L K

其中，h_i是有序的特征直方图值序列H中的第i个值，1≤i≤μ，μ表示特征序列

中取不同值的元素个数，Ω_j是第j个聚类集合，1≤j≤K，K表示聚类个数，

是第θ次迭代时第j个聚类集合的聚类中心，是第θ次迭代时第l个聚类集合的聚类中心，1≤l≤K；

D5)更新K个聚类集合的聚类中心：

$D_j^{(\mathrm{\&theta;}+1)}=\underset{h_i\&Element;{\mathrm{\&Omega;}}_j}{\mathrm{\&Sigma;}}{h}_i/\mathrm{num}\left({\mathrm{\&Omega;}}_j\right),$ 1≤j≤K

其中，

是第θ+1次迭代时第j个聚类集合Ω_j的聚类中心，1≤j≤K，K表示聚类个数，h_i是有序的特征直方图值序列H中的第i个值，1≤i≤μ，μ表示特征序列

中取不同值的元素个数，num(Ω_j)表示第θ次迭代时第j个聚类集合Ω_j中的元素个数；然后，将迭代次数θ加1；

D6)重复步骤D4)～D5)直到找到满足条件

的聚类集合Ω₁，L，Ω_K或者使迭代次数θ达到100次，其中，h_i是有序的特征直方图值序列H中的第i个值，1≤i≤μ，μ表示特征序列

中取不同值的元素个数，Ω_j是第j个聚类集合，

是第θ+1次迭代时第j个聚类集合的聚类中心，1≤j≤K，K表示聚类个数，将聚类集合Ω₁，L，Ω_K组成一个集合，得到分类序列Ω＝{Ω₁，L，Ω_K}，Ω_K是分类序列Ω中的第K个值。

步骤E：待检测嵌入块序列E^S的特征值与分类序列Ω进行比较，提取水印序列R，并用嵌入级别数L修改待检测嵌入块序列E^S的特征值，获得重构图像I^R。

E1)根据分类序列Ω＝{Ω₁，L，Ω_K}，确定K-1个类边界点Q₁，L，Q_K-1：当K＝3时，类边界点表示为：Q₁＝max{Ω₁}+1，Q₂＝min{Ω₃}-1；当K＝2时，类边界点表示为Q₁＝max{Ω₁}+1，其中，Ω₁为分类序列Ω中第1个聚类集合，Ω₃为分类序列Ω中第3个聚类集合，max{Ω₁}是取第1个聚类集合Ω₁中的最大值，min{Ω₃}是取第3个聚类集合Ω₃中的最小值；

E2)将待检测嵌入块序列E^S中第i个块的特征值

与K-1个类边界点Q₁，L，Q_K-1进行比较，提取第i位水印r_i：

E2a)当K＝3时，提取规则为：

其中，是待检测嵌入块序列中第i个块的特征值，Q₁是分类序列Ω中第1个类边界点，Q₂是分类序列Ω中第2个类边界点；

E2b)当K＝2时，如果F_l＝0且F_r≠0，提取规则为：

其中，F_l是峰值对W的左频率，F_r是峰值对W的右频率，

是待检测嵌入块序列中E^S第i个块的特征值，Q₁是分类序列Ω中第1个类边界点；

E2c)当K＝2时，如果F_l≠0且F_r＝0，提取规则为：

E3)重复步骤E2)直到得到提取的水印序列R＝{r₁，Kr_i，Kr_f′}，其中，r_i是水印序列R中的第i位水印，r_i∈{0，1}，1≤i≤f′，f′表示水印序列R的长度；

E4)用嵌入级别数L修改待检测嵌入块序列E^S中第i个块的特征值，获得重构的嵌入块：

E4a)当K＝3时，修改规则为：

其中，ξ′_i是重构的嵌入块序列中第i个块的特征值，

是待检测嵌入块序列中第i个块的特征值，L是嵌入级别数，Q₁是分类序列Ω中第1个类边界点，Q₂是分类序列Ω中第2个类边界点；

E4b)当K＝2时，如果F_l＝0且F_r≠0，修改规则为：

其中，F_l是峰值对W的左频率，F_r是峰值对W的右频率，ξ′_i是重构的嵌入块序列中第i个块的特征值，是待检测嵌入块序列中第i个块的特征值，L是嵌入级别数，Q₁是分类序列Ω中第1个类边界点；

如果F_l≠0且F_r＝0，修改规则为：

E4c)根据步骤E4a)～E4b)中修改第i个块的特征值的不同规则，修改待检测嵌入块序列中第i个块的像素值：

E4c1)当

时，将待检测嵌入块序列E^S中第i个块中所有像素的灰度值保持不变；

E4c2)当

时，将待检测嵌入块序列E^S中第i个块的两个像素集合

和

中像素的灰度值修改为：

式中，

是重构后第i个块的像素集合

中第k个像素的灰度值，

是待检测嵌入块序列E^S中第i个块的像素集合

中第k个像素的灰度值，L是嵌入级别数，是重构后第i个块的像素集合

中第k个像素的灰度值，

是待检测嵌入块序列E^S中第i个块的像素集合

中第k个像素的灰度值，t为表示待检测图像时所用的二进制位数，取8；

E4c3)当

时，将待检测嵌入块序列E^S中第i个块的两个像素集合

和

中像素的灰度值修改为：

E5)重复步骤E4)直到将待检测嵌入块序列E^S中的所有块都重构出来，并根据嵌入块序列E在原始图像块序列B中的位置信息S，将重构的块和未重构的块进行组合，得到重构图像I^E。

通过上述步骤A～步骤E可以实现水印的提取，从待检测的图像I^S中提取出二值的伪随机水印序列R，并且得到重构图像I^E。

本发明的优点可通过以下仿真实验进一步说明：

将本发明方法分别与现有基于直方图旋转的方法和现有基于直方图分布约束的方法进行实验对比，实验中首先用这些方法分别将二值的伪随机水印序列嵌入到大小为512×512×8的灰度图像中得到含水印图像，然后对生成的含水印图像进行JPEG和JPEG2000的有损压缩得到降质的含水印图像，最后再利用这些方法从降质的含水印图像中提取水印进行容量、视觉质量和鲁棒性的实验。

1.实验条件与实验说明

实现本发明的软件环境为美国Mathworks公司开发的MATLAB 2009a，实验中大量的灰度图像包括来自CVG-UGR图像数据库的自然图像和来自DICOM样本图像集与OsiriX网站的医学图像。本发明的一些记号为：将基于直方图旋转的方法记为HR，将基于直方图分布约束的方法记为HDC，将本发明方法记为HN。实验中三种方法的原始图像和待检测图像的分块大小均为8×8，并且嵌入强度相同，即HR中嵌入级别数为8，HDC中阈值为8，HN方法中嵌入级别数为16，而且在HDC中采用BCH(15，11，1)作为误差纠错码。

2.实验内容

实验1：容量实验

本发明以客观指标位/像素bpp作为评判依据，测试现有HR、HDC和本发明HN三种方法的容量，bpp＝Q/PX，其中Q是指在原始图像中最多能嵌入的水印位数，PX是指在原始图像中像素的个数，在本实验中PX＝512×512＝262144。在原始图像大小确定的情况下，bpp越大，表示嵌入水印的位数越多，即容量越大。本实验在自然图像数据库上进行测试，实验结果如表1所示。

表1.不同方法的容量值

由表1结果可见，由于没有采用误差纠错码，本发明方法的容量明显高于HDC，克服了HDC方法容量低的缺陷，同时，通过块大小和阈值控制容量的方式使得本发明方法的容量与HR相当。

实验2：视觉质量实验

本发明以客观指标峰值信噪比PSNR作为评判依据，测试三种方法在原始图像以最大容量嵌入情况下含水印图像的视觉质量，其中PSNR表示为

$\mathrm{PSNR}=10\log \left(\frac{M\&times;N\&times;{255}^2}{\underset{i=1}{\overset{M}{\mathrm{\&Sigma;}}}\underset{j=1}{\overset{N}{\mathrm{\&Sigma;}}}{(I(i,j)-I^E(i,j))}^2}\right),$

式中，M×N是原始图像大小，I(i，j)是原始图像在第i行第j列的像素值，I^E(i，j)是含水印图像在第i行第j列的像素值。本实验在自然图像和医学图像数据库上进行测试，实验结果如图3和表2所示，其中图3(a)为HR方法生成的含水印图像，图3(b)为本发明方法生成的含水印图像。

表2.不同方法的PSNR(dB)值

从表2的实验结果可以看出，本发明的PSNR明显高于HR和HDC方法，特别是在医学图像数据库上，本发明方法的PSNR比HR有了显著的提高，其中PSNR越高，代表含水印图像的视觉质量越好，说明本发明的有效性。

由图3结果可见，本发明方法由于采用块筛选的方式有效解决了水印嵌入过程中的像素溢出问题，从而去除了HR方法含水印图像中的“椒盐”噪声，克服了HR方法视觉质量差的缺点。

实验3：鲁棒性实验

本发明进行鲁棒性实验的具体过程是：在水印嵌入过程中，将长度等于固定值的水印序列重复嵌入到原始图像中，重复的次数等于

其中V是重复嵌入的次数，

是向下取整运算，即取小于或等于Q/TR的最大整数，Q是指在原始图像中最多能嵌入的水印位数，TR是固定值，实验中取100；在水印提取过程中，将提取出来的水印序列分成长度等于TR的V个子序列，再将这些子序列按位分别组成一组，将每组中的V个值采用多数表决的方式进行判定，最后得到提取的水印序列。

本发明在JPEG和JPEG2000两种常见的有损压缩攻击下测试三种方法的鲁棒性，实验中JPEG压缩的质量因子取9个值：20，30，40，50，60，70，80，90，100，分别在每个质量因子下测试三种方法的鲁棒性；JPEG2000压缩的压缩率取10个值：0.2，0.4，0.6，0.8，1.0，1.2，1.4，1.6，1.8，2.0，分别在每个压缩率下测试三种方法的鲁棒性，其中JPEG压缩的质量因子和JPEG2000压缩的压缩率越小，压缩强度越大。

本发明以客观标准存活率SL作为评判依据，测试三种方法的鲁棒性，SL表示为SL＝|φ-δ|/φ，其中0≤SL＜1，φ是最弱压缩强度，JPEG压缩时φ取100，JPEG2000压缩时φ取2.0，δ是恰可提取压缩强度，即当压缩强度小于或者等于δ时，嵌入的水印序列恰能被正确提取出来，SL越大，鲁棒性越强，反之亦然。同时，本发明以客观标准鲁棒置信度RCR作为评判依据，测试三种方法鲁棒性的可靠性，RCR表示为RCR＝Φ/Θ，其中Φ是在一个数据库上使用最弱压缩强度能正确提取水印序列的原始图像个数，Θ是该数据库中全部原始图像的个数，RCR越大，鲁棒性的可靠性越高，反之亦然。

在医学图像数据库上，可靠性的实验结果如表3所示，其中在JPEG压缩下，本发明方法与HR方法的鲁棒性均为0.5，在JPEG2000压缩下，HR方法的鲁棒性是0.4，而本发明方法的鲁棒性是0.9。从表3的实验结果可见，在相同或者高于HR方法鲁棒性的情况下，本发明方法的鲁棒置信度要优于HR方法，解决了HR方法鲁棒性不可靠的问题。

表3.不同方法的鲁棒置信度

在自然图像数据库上，可靠性的实验结果如表4所示，其中在JPEG压缩下，本发明方法与HR方法的鲁棒性均为0.5，在JPEG2000压缩下，两种方法的鲁棒性是0.7。由表4结果可见，在相同的鲁棒性下，本发明方法的鲁棒置信度要高于HDC方法。

表4.不同方法的鲁棒置信度

综上，本发明克服了HR方法含水印图像有“椒盐”噪声、视觉质量差和HDC方法容量低的缺点，解决了它们鲁棒性不可靠的问题，提高了鲁棒可逆图像水印方法的综合性能。

Claims (9)

1.一种基于直方图邻域的鲁棒可逆水印嵌入方法，包括如下步骤：

(1)将原始图像I分成大小相同且互不重叠的块，丢弃不能被完全划分的块，得到一组原始图像块序列B；

(2)将原始图像块序列B进行筛选，得到候选块序列A；

(3)计算候选块序列A的特征值并生成特征直方图，在特征直方图中确定峰值对W；

(4)用峰值对W在候选块序列A中查找嵌入块，得到一组嵌入块序列E以及嵌入块序列E在原始图像块序列B中的位置信息S；

(5)确定嵌入级别数L，L为正整数，并利用伪随机数产生函数生成二值的伪随机水印序列P；

(6)用嵌入级别数L修改嵌入块序列E的特征值，将水印序列P嵌入到嵌入块序列E中，得到含水印图像I^E：

(6a)从水印序列P中选择第j位水印p_j，p_j∈{0，1}，1≤j≤f，f表示水印序列P的长度；

(6b)修改嵌入块序列E中第i个块的特征值ξ_i，将水印p_j进行嵌入，嵌入规则为：

式中，

是嵌入水印后嵌入块序列E中第i个块的特征值，ξ_i是嵌入水印前嵌入块序列E中第i个块的特征值，1≤i≤σ，σ为嵌入块序列E中块的个数，p_j是水印序列中的第j位水印，L是嵌入级别数，ξ_l是峰值对W的左值，ξ_r是峰值对W的右值；

(6c)重复步骤(6a)～(6b)直到将水印序列P中的所有位水印都嵌入到原始图像I中，得到含水印图像I^E。

2.根据权利要求1所述的鲁棒可逆水印嵌入方法，其中步骤(2)所述的将原始图像块序列B进行筛选，按如下步骤进行：

(2a)将原始图像块序列B中第i个块的所有像素分成两个大小相同的像素集合

和

1≤i≤v，v为原始图像块序列B中块的个数：

$Z_a^i=\left\{Z^{(x,y)}\right|Z^{(x,y)}\&Element;Z,\mathrm{mod}(x+y,2)=0\}$

$Z_b^i=\left\{Z^{(x,y)}\right|Z^{(x,y)}\&Element;Z,\mathrm{mod}(x+y,2)\&NotEqual;0\},$

式中，Z^(x，y)是第i个块中位于第x行、第y列的像素的灰度值，1≤x≤m，1≤y≤n，m×n是块大小，Z是第i个块中所有像素灰度值的集合，mod(·)是取余数运算，mod(x+y，2)是计算x+y除以2所得到的余数；

(2b)从得到的两个像素集合

和中分别按顺序选取一个像素，用这两个像素的灰度值构成一个像素对

其中

代表从

中所选择的第k个像素的灰度值，

代表从

中所选择的第k个像素的灰度值，1≤k≤m×n/2；

(2c)按如下方式对第i个块中的所有像素对

1≤k≤m×n/2进行判定：

$\left\{\begin{array}{c}{Z}_a^{(i,k)}<L,\\ Z_b^{(i,k)}>2^t-1-L\end{array}\right.$ 或者 $\left\{\begin{array}{c}{Z}_a^{(i,k)}>2^t-1-L,\\ Z_b^{(i,k)}<L\end{array}\right.,$

式中，L是嵌入级别数，t为表示原始图像时所用的二进制位数，取8；若第i个块中的所有像素对

1≤k≤m×n/2都不满足上述条件，认为第i个块为候选块；

(2d)重复步骤(2a)～(2c)直到将原始图像块序列B中每个块都判断完为止，得到一组候选块序列A。

3.根据权利要求1所述的鲁棒可逆水印嵌入方法，其中步骤(3)所述的计算候选块序列A的特征值并生成特征直方图，在特征直方图中确定峰值对W，按如下步骤进行：

(3a)计算候选块序列A中每个块的特征值，组成特征序列Ψ＝(ξ₁，Kξ_i，Kξ_ω)，其中，ξ_i是第i个块的特征值，1≤i≤ω，ω为候选块序列A中块的个数；

(3b)根据特征序列Ψ中的取值，生成候选块序列A的特征直方图，该直方图包括值序列X＝{x₁，L x_i，L x_η}和频率序列F＝{F₁，L F_i，L F_η}，其中，x_i是值序列X中的第i个值，F_i是频率序列F中的第i个值，1≤i≤η，η表示特征序列Ψ中取不同值的元素个数；

(3c)在特征直方图的频率序列F中找出频率取最大值和取次最大时所对应的值序列X中的特征值，生成峰值对W＝[(ξ_l，F_l)，(ξ_r，F_r)]，其中，ξ_l是峰值对W的左值，ξ_r是峰值对W的右值，F_l是峰值对W的左频率，F_r是峰值对W的右频率。

4.根据权利要求1所述的鲁棒可逆水印嵌入方法，其中步骤(4)所述的用峰值对W在候选块序列A中查找嵌入块，按如下步骤进行：

(4a)计算候选块序列A中第i个块的特征值ξ_i与峰值对W的左值ξ_l之间的距离d(ξ_i，ξ_l)＝|ξ_i-ξ_l|，以及特征值ξ_i与峰值对W的右值ξ_r之间的距离d(ξ_i，ξ_r)＝|ξ_i-ξ_r|，1≤i≤ω，ω为候选块序列A中块的个数；

(4b)当ξ_i≤ξ_l且d(ξ_i，ξ_l)小于或等于给定的阈值T或者当ξ_i≥ξ_r且d(ξ_i，ξ_r)小于或等于给定的阈值T时，认为第i块是嵌入块，记录该块在原始图像块序列B中的位置；

(4c)重复步骤(4a)～(4b)，直到候选块序列A中所有块都判断完为止，得到一组嵌入块序列E以及嵌入块序列E在原始图像块序列B中的位置信息S。

5.一种基于直方图邻域的鲁棒可逆水印提取方法，包括如下步骤：

1)将待检测图像I^S分成大小相同且互不重叠的块，得到一组待检测图像块序列C；

2)对待检测图像块序列C进行筛选，得到一组待检测嵌入块序列E^S；

3)计算待检测嵌入块序列E^S的特征值并生成特征直方图，将特征直方图中的值序列进行排序，得到一组有序的特征直方图值序列H；

4)利用K-means聚类算法将有序的特征直方图值序列H进行划分，得到一组分类序列Ω：

4a)利用水印嵌入过程中得到的峰值对W的左频率F_l和右频率F_r的取值，对K-means聚类算法进行初始化，初始化的值包括聚类个数K、聚类集合Ω₁，L，Ω_K，聚类中心D₁，L，D_K和迭代次数θ；

4b)按照特征直方图值序列H中每个值与K个聚类中心D₁，L，D_K之间的距离，将每个值分配到最相近的一个聚类集合中；

4c)更新K个聚类中心和迭代次数θ；

4d)重复步骤4b)～4c)直到找到满足条件

的聚类集合Ω₁，L，Ω_K或者使迭代次数θ达到100次，其中，h_i是特征直方图值序列H中的第i个值，1≤i≤μ，μ表示特征序列

中取不同值的元素个数，Ω_j是第j个聚类集合，是第θ+1次迭代时第j个聚类集合的聚类中心，1≤j≤K，K表示聚类个数，将聚类集合Ω₁，L，Ω_K组成一个集合，得到分类序列Ω＝{Ω₁，L，Ω_K}，Ω_K是分类序列Ω中的第K个值；

5)将待检测嵌入块序列E^S的特征值与分类序列Ω进行比较，提取水印序列R，并用嵌入级别数L修改待检测嵌入块序列E^S的特征值，获得重构图像I^R。

6.根据权利要求5所述的鲁棒可逆水印提取方法，其中步骤2)所述的对待检测图像块序列C进行筛选，得到一组待检测嵌入块序列E^S，按如下步骤进行：

2a)利用水印嵌入过程中得到的嵌入块序列在原始图像块序列中的位置信息S＝(s₁，Ks_i，Ks_v)，s_i是向量S中的第i个值，对待检测图像块序列C中的第i个块进行判定：如果s_i＝1，认为该块是待检测嵌入块，如果s_i＝0，认为该块不是待检测嵌入块，其中，1≤i≤v，v为原始图像块序列中块的个数，该个数与待检测图像块序列C中的块个数相同；

2b)重复步骤2a)直到待检测图像块序列C中的所有块都判断完为止，得到一组待检测嵌入块序列E^S。

7.根据权利要求5所述的鲁棒可逆水印提取方法，其中步骤3)所述的计算待检测嵌入块序列E^S的特征值并生成特征直方图，将特征直方图中的值序列进行排序，按如下步骤进行：

3a)计算待检测嵌入块序列E^S中每个块的特征值，组成特征序列

其中，

是第i个块的特征值，1≤i≤σ′，σ′为待检测嵌入块序列E^S中块的个数；

3b)根据特征序列中的取值，生成待检测嵌入块序列E^S的特征直方图，该直方图包括值序列

和频率序列

其中，

是值序列

中的第i个值，是频率序列

中的第i个值，1≤i≤μ，μ表示特征序列

中取不同值的元素个数；

3c)将值序列

按照取值进行递增排序，得到一组有序的特征直方图值序列H＝[h₁，L h_i，L h_μ]，其中，h_i是有序的特征直方图值序列H中的第i个值，1≤i≤μ，μ表示特征序列

中取不同值的元素个数，h₁＜L＜h_i＜L＜h_μ。

8.根据权利要求5所述的鲁棒可逆水印提取方法，其中步骤5)所述的将待检测嵌入块序列E^S的特征值与分类序列Ω进行比较，提取水印序列R，按如下步骤进行：

5a)根据分类序列Ω＝{Ω₁，L，Ω_K}，确定K-1个类边界点Q₁，L，Q_K-1：当K＝3时，类边界点表示为：Q₁＝max{Ω₁}+1，Q₂＝min{Ω₃}-1；当K＝2时，类边界点表示为Q₁＝max{Ω₁}+1，其中，Ω₁为分类序列Ω中第1个聚类集合，Ω₃为分类序列Ω中第3个聚类集合，max{Ω₁}是取第1个聚类集合Ω₁中的最大值，min{Ω₃}是取第3个聚类集合Ω₃中的最小值；

5b)将待检测嵌入块序列E^S中第i个块的特征值

与K-1个类边界点Q₁，L，Q_K-1进行比较，提取第i位水印r_i：

5b1)当K＝3时，提取规则为：

其中，

是待检测嵌入块序列E^S中第i个块的特征值，Q₁是分类序列Ω中第1个类边界点，Q₂是分类序列Ω中第2个类边界点；

5b2)当K＝2时，如果F_l＝0且F_r≠0，提取规则为：

其中，F_l是峰值对W的左频率，F_r是峰值对W的右频率，

是待检测嵌入块序列中E^S第i个块的特征值，Q₁是分类序列Ω中第1个类边界点；

5b3)当K＝2时，如果F_l≠0且F_r＝0，提取规则为：

5c)重复步骤5b)直到得到提取的水印序列R＝{r₁，Kr_i，Kr_f′}，其中，r_i是水印序列R中的第i位水印，1≤i≤f′，f′表示水印序列R的长度。

9.根据权利要求5所述的鲁棒可逆水印提取方法，其中步骤5)所述的用嵌入级别数L修改待检测嵌入块序列E^S的特征值，获得重构图像I^R，按如下步骤进行：

5A)用嵌入级别数L修改待检测嵌入块序列E^S中第i个块的特征值，获得重构的嵌入块：

5A1)当K＝3时，修改规则为：

其中，ξ′_i是重构的嵌入块序列中第i个块的特征值，

是待检测嵌入块序列中第i个块的特征值，L是嵌入级别数，Q₁是分类序列Ω中第1个类边界点，Q₂是分类序列Ω中第2个类边界点；

5A2)当K＝2时，如果F_l＝0且F_r≠0，修改规则为：

其中，F_l是峰值对W的左频率，F_r是峰值对W的右频率，ξ′_i是重构的嵌入块序列中第i个块的特征值，

是待检测嵌入块序列中第i个块的特征值，L是嵌入级别数，Q₁是分类序列Ω中第1个类边界点；

当K＝2时，如果F_l≠0且F_r＝0，修改规则为：

5A3)根据步骤5A1)～5A2)中修改第i个块的特征值

的不同规则，修改待检测嵌入块序列中第i个块的像素值；

5B)重复步骤5A)直到将待检测嵌入块序列E^S中的所有块都重构出来，并根据嵌入块序列E在原始图像块序列B中的位置信息S，将重构的块和未重构的块进行组合，得到重构图像I^E。

Images