基于K-means聚类和Harris特征的水印算法∗

郑秋梅 杨亚男 曹宝琴 王风华 张萌萌

（中国石油大学（华东）计算机与通信工程学院 青岛 266580）

1 引言

随着互联网、多媒体技术的迅速发展，使得多媒体内容分发、通讯、再现变得容易，随之而来也出现了很多问题，如内容的篡改、非法复制、发行、再传输。数字水印技术［1～2］是解决这些问题的有效手段之一。数字水印是以不易察觉的方式将水印嵌入作品所有权版权信息（企业标识、作品信息等），来保护版权所有者的合法利益。目前数字水印算法大多数以灰度图像为研究载体，而彩色图像在实际的应用场景中占主导地位。目前彩色图像数字水印算法大多利用人类视觉系统的特性，在载体图像的频域系数中嵌入水印，然而由于频域系数没有几何不变性，不能有效抵抗几何攻击，因此，设计彩色图像能够抵抗几何攻击的水印算法成为研究热点和难点。

近几年抗几何攻击的水印算法主要采用第二代数字水印技术［3～4］。第二代数字水印技术利用图像内在特征进行水印嵌入在抗几何攻击中具有良好的鲁棒性。其基本思想是利用提取图像中具有几何不变性的特征点，利用这些特征点来构造特征区域，将水印信息嵌入其中，从而提高水印算法对几何攻击的鲁棒性。目前针对抗几何攻击，国内外学者利用Harris 角点具有RST 不变性的特性做了大量研究。Qi等［5］提出，图像内容由基于图像纹理的自适应Harris 角检测器获得的重要特征点表示，这些重要的特征点具有几何意义，因此能够利用基于Delaunay三角形匹配方法去同步攻击。陈青等［6］利用Harris 算子提取并筛选稳定的特征点，在特征区域的小波域嵌入水印，该算法抗旋转攻击鲁棒性较好，但水印图像和载体图像像素大小受特征区域个数限制。周广州等［7］利用Harris算子提取出归一化图像的特征点，选取部分稳定特征点确定特征区域，对特征区域进行一次小波分解，利用分块加权奇异值分解法提高水印的鲁棒性，该算法对剪切攻击鲁棒性良好，但对中值滤波、缩放攻击、旋转攻击等水印攻击鲁棒性不佳。

单一地使用Harris 算子提取特征点，通过这些特征点选定的特征区域，特征点个数较少，当图像旋转角度较大时，遭受攻击后点位移或剪切时，特征区域不稳定，水印无法恢复。本文根据人眼视觉特性，采用彩色图像作为载体图像，利用具有旋转不变性的Harris 算子，提出一种新的基于聚类算法的特征构造算法，使得到的特征区域不但尽可能多地包含稳定的特征点又互不重叠，采用DC 系数自适应的水印嵌入算法，保证了对常规图像处理攻击的鲁棒性。

2 背景及相关理论知识

2.1 Harris角点检测

Harris 角点检测函数计算简单，同时该响应函数具有旋转不变性。Harris［8～9］在1988 年改良并优化Moravec 算子，考虑角点与方向的差值，提出了Harris 角点检测算法。Moravec 算法［10］的窗口移动方向都是45°或45°的整数倍，而Harris 算法将窗口沿各个方向进行移动，计算出各个方向上的像素梯度变化，定义灰度值的变化如式（1）：

其中，Iu，v为图像的灰度函数，wu，v=e-(u2+v2)/ο2为窗口函数。

对于局部微小的移动量[x，y] ，可以近似得到下面的表达：

其中，k是一个经验常量，一般k∈[0 .04，0.06] ，Tr(M)和Det(M)分别是矩阵M的行列式和迹。

R值与M的特征值有关，当R为大数值正数、大数值负数和小数值，对应检测到的是角点、边缘和平坦区。图1 为对Lena 图和Baboon 图进行Har⁃ris特征检测后得到的结果。

图1 Harris特征检测结果

由图1 可以看出，Harris 角点的检测算法提取出的角点数目较多且密集，分布较广，聚类算法构造特征区域具有局限性，因此本文采用一种加权Harris 角点检测算法［11］，给每个检测到的角点分配权值，通过权值大小来进行角点筛选。

将式（5）中计算得到的R值归一化，作为其初始权值，Harris 角点检测算法提取出的角点作为窗口的中心点，每一个角点用窗口（尺寸为7×7）进行遍历。根据式（6）计算中心角点的权值。

其中value 为权值，(x，y)为中心角点的坐标，(x'，y')为窗口覆盖区域坐标。

根据value 的大小对提取到的角点进行划分，权值大的角点比权值小的角点更为重要，利用特征权值对特征点进行初步筛选后，采用K-means聚类算法构造特征区域。

2.2 K-means聚类算法

K-means 是解决众所周知的聚类问题最简单的无监督学习算法之一。其基本思想［12］为将n个数据点划分为K个不相交的类，其类簇中心为μi，把数据点划分至与μi距离最小的类簇中，以便最小化等式（7）的平方误差函数：

K-means算法详细流程如下。

步骤1初始化质心μi，i=1，2，...k；

步骤2按照距离最近原则，将数据点划分到最近的聚类中心那一类中；

步骤3当分配完所有点后，计算K 个质心的位置；

步骤4循环步骤2 和步骤3，直至质心不再改变。

对K 值和初始质心的选取方面，K-means 算法需要选取K 个初始质心，常见的方法是随机选取，但这种方法得到的聚类划分效果较差；K-means算法选取初始质心需要不停地计算距离平方和，以得到最佳聚类中心，若数据集里面的数据量很大时，需要耗费大量的时间。本文采用与层次聚类结合的方法［13］，首先对特征点进行层次聚类，根据层次聚类结果选取K值和初始质心。

3 基于K-means 聚类和Harris 角点的特征区域构造算法

不同特征区域内可能会含有相同的特征点，若构造方法不得当不仅会导致特征区域重叠反复，还容易导致水印信息的交错相互影响，因此利用K-means聚类算法设计特征区域构造算法，使得到的特征区域不但包含稳定的特征点又互不重叠。

本文提出一种新的基于K-means 聚类和Har⁃ris角点的特征区域构造算法，通过加权Harris角点检测提取图像的特征点，选取响应值R较大的特征点，利用K-means聚类算法设计特征区域构造算法，使得到的特征区域不但尽可能多的包含稳定的特征点又互不重叠。

构造算法具体流程如下。

步骤1通过加权Harris 角点检测算法给每个检测到的角点分配权值value，权值较大的点周围的角点分布密集，权值较小的点周围的角点数量少；

步骤2根据需要构造的特征区域数TN 选取合适的K 值，根据角点value 值的大小初步筛选后，利用层次聚类法对特征点初次聚类，得到K 个层次聚类中心，以求得的K 个中心点作为K-means聚类初始中心，进行K-means 聚类得到最终的K 个K-means聚类中心；

步骤3根据K 个聚类中心与嵌入区域半径大小r 对中心角点进行筛选，去掉区域超过载体图片大小的角点，得到K1个中心点；

步骤4比较距离中心点最近的角点的value值大小，value 值最大中心点作为区域中心进行特征区域构造，以此方法构造特征区域内焦点数目，权值较大的点周围的角点分布密集；

步骤5删除位于步骤4 构造的特征区域内的中心点，得到新的中心点序列；

步骤6判断已构造的特征区域数是否为TN，若小于则执行步骤4，若等于，则构造完成，使得到的特征区域不但尽可能多的包含稳定的特征点又互不重叠。

4 水印嵌入和提取算法

4.1 水印嵌入算法

水印的嵌入流程图如图2，本文嵌入算法基本流程如下。

步骤1水印预处理

大小为m×m的二值图像作为水印图像，采用Arnold 变换和Logistic 混沌序列［14］相结合的方法对水印进行预处理，将Arnold 变换次数n、Logistic 映射的迭代初始值x0、迭代控制参数μ以及去除序列项数t作为密钥保存。

今天，三北的生态状况依然脆弱，三北的生态建设任务依然任重道远。也许，就在未来的某一天，沙暴或者沙魔还会来袭，如果防沙不力，已有的成果可能葬身沙海，被无情地埋葬。

步骤2提取Harris角点构造局部特征区域

构造嵌入区域，将构造特征区域的中心点坐标信息、嵌入区域半径r 及载体图像尺寸S 作为密钥保存，此外选取权值大的u 个特征点信息作为密钥保存。然后将求得的局部特征区域进行合并，将合并后的区域作为水印的嵌入区域。

图2 本文算法水印嵌入流程图

步骤3水印位判断取反

将上述得到的嵌入区域B 分量作为水印嵌入区域，下采样为EB1 和EB2 两部分。分别对区域EB1 和区域EB2 进行8×8 分块DCT 变换，得到DC系数分别记为DCi与DCj，根据式（8）对要嵌入的1bit水印信息W进行取反判断：

其中，α为水印的嵌入强度，将纹理不同的8×8子块经过DCT变换得到两个AC系数带入式（10），判断子块类型，然后根据式（11）选取对应的嵌入强度。这里阈值β的选取参照了Wang［15］中的经验阈值0.05。具体判断方法如下所示：

式中的λ1，λ2取经验阈值：λ1=50，λ2=65。

其中，αsmooth＜αedge＜αtexture，通过实验验证，对于RGB 模 型 设 置 为：αsmooth=0.07 ，αtexture=0.2 ，αedge=0.1。

嵌入区域图像的RGB 模型的G 分量进行和B分量相同的操作，G 分量嵌入另一半水印，然后进行DCT逆变换，得到最终嵌入水印的图像。

4.2 水印提取算法

本文在水印提取前，首先利用密钥中保存的权值较大的特征点信息进行几何校正［16］，利用密钥中保存的特征区域中心点坐标构造特征区域，对校正后的图像提取水印，其具体步骤如下。

步骤1通过加权Harris 角点检测算法提取遭受攻击后的图像特征，保存权值较大的u 个特征点信息，与密钥中保存的特征进行匹配；

步骤2图像几何校正，对图像缩放攻击的矫正直接根据密钥中保存的载体图像尺寸S 进行缩放攻击校正。对图像旋转攻击的矫正通过式（12）计算逆旋转角度Θ 进行旋转攻击校正。计算公式如下：

θ为旋转角度。

步骤3根据密钥中保存的区域中心点构造局部区域，转化为RGB 模型，提取B 分量下采样记为TB1和TB2；

步骤4分别对TB1与TB2进行8×8分块DCT变换，得到对应的DC 系数DCi与DCj，根据式（13）进行水印提取。

根据密钥中保存的取反坐标序列，对提取到的水印信息进行取反处理，得到一半水印信息。

步骤5取RGB 模型中的G 分量，下采样为TG1 和TG2，分别对TG1 与TG2 进行8×8 分块DCT变换，得到对应的DC 系数DCj与DCj，利用式（13），根据密钥中保存的取反坐标序列，对提取到的水印信息进行取反处理，提取到另一半水印。

步骤6将步骤4和步骤5水印信息合并，利用密钥进行解密，合并得到最终的水印图像。

5 实验结果分析

实验使用Matlab 2017b 作为验证平台，实验中采用512 pixel×512 pixel 彩色Lena 图像作为载体图像，选择30 pixel ×30 pixel的二值图像作为嵌入的水印图像，图3为实验使用的载体图像及水印图像，图4为原始图像及嵌入水印后的图像。

采用峰值信噪比（PSNR）作为算法不可见性的评价指标［17］，通过实验仿真，本文算法的PSNR 值为45dB，远远高于人眼可感知程度的33dB，具有很好的不可见性。

图3 载体图像和水印图像

图4 原始图像和嵌入水印后的图像

1）用归一化系数（NC）值客观评价原始图像与提取水印后图像的相关度［18］。常规信号攻击包括不同强度的压缩、噪声和滤波等，本文算法和马婷等［19］提取出的水印信息与原始水印信息的NC值见表1，实验数据表明本文提出的水印算法对常规信号攻击的鲁棒性上要明显优于对比算法。

表1 常规图像攻击下的NC值

2）旋转攻击。本文算法和陈青等［20］提取出的水印信息与原始水印信息的NC 值见表2，本文算法经旋转攻击后提取的NC 值接近1，实验数据表明本文提出的水印算法在抗旋转攻击下的鲁棒性上要明显优于对比算法。

表2 旋转攻击下的NC值

3）缩放攻击。采用双三次插值方法对水印进行同等比例缩放攻击（缩放比例分别为0.5，0.75，1.25，1.5 倍）。表3 为本文算法和文献［20］提取出的水印信息与原始水印信息的NC 值，对比文献［20］可知，本文算法在抵抗缩放攻击更有优势。

表3 缩放攻击下的NC值

6 结语

本文针对彩色图像数字水印算法对几何攻击鲁棒性较差的问题，提出一种基于K-means聚类和Harris 特征的第二代数字水印算法。选取具有旋转不变性的Harris 角点来构造水印嵌入区域。通过对特征的区域构造算法进行改进，采用K-means聚类算法确定嵌入区域的中心，保证了嵌入区域内特征点的最大化。由Matlab仿真实验分析，本文算法的PSNR 值达到45dB，远远高于人眼可感知程度的33dB，具有很好的不可见性；此外本文算法在不同图像攻击下，计算得到的NC 值接近1，不仅对常规图像处理具有鲁棒性，能够很好地抵抗几何攻击。而相较于大部分的基于特征区域的水印算法，本文提出的算法具有更高的水印容量。