Cover-Lossless Robust Watermarking【文章学习】

概要----本文主要针对“cover-lossless robust watermarking”的主要原理理解，代码解读，特别关注Zernike距的理解，没有结果部分。

目录

摘要

相关前人工作

准备工作

几何攻击

泽尼克Zernike矩

前人基于Zernike矩的方法

无损鲁棒水印方法（附代码）

A.新的水印嵌入策略（水印嵌入到Zernike距中）

1.先计算Zernike距

2.将比特嵌入选定Zernike距

3.可逆Zernike变换

4.舍入误差的计算（第三个误差）

5.补偿信息的可逆嵌入

B.完整性认证

C.提取水印和恢复图像（前提是不受攻击）

D.在攻击在的水印提取

---

摘要

无损鲁棒水印（cover-lossless robust watermarking）---在无攻击的条件下恢复原始的图像（再也不用担心加水印而丢失图像了！）

前人所提出的方法都为加性方法：JPEG compression, low-pass filtering and Gaussian additive noise。。。不可以抗几何变形攻击（geometric deformations：rotation and scaling）。The main reason is that those exploited robust features are related to the pixel position.

这篇文章的方法：
将水印有效地嵌入到低阶泽尼克矩（low-order Zernike moments）中，并将引起的失真可逆隐藏，作为原始图像恢复的补偿信息。

鲁棒特征：泽尼克矩（low-order Zernike moments）不仅对旋转不变性，而且抗缩放（尺度不变性）和平移不变性。
.....

实验结果表明，所提出的鲁棒水印方法可以有效地减少补偿信息，并且对缩放、旋转、JPEG压缩等类噪声操作具有较强的鲁棒性。在无攻击下，可以无损恢复图像。

相关前人工作

分大概两类 可逆水印和可逆鲁棒水印

可逆水印（reversible watermarking）：利用像素间的冗余来嵌入和提取水印，不损伤像素（可逆=无损）。应用：脆弱身份验证；限制：对内容保存的敏感性？

可逆鲁棒水印（又称 原图无损鲁棒水印）：
基于直方图位移：通过在空间域修改图像的差分直方图来实现可逆数据隐藏。抗加性噪声通过分块，并通过移动直方图两三个区域来嵌入水印，其中鲁棒性通过区域距离获得，但是水印的直方图与块（blocks）的位置有关，因此对于几何变形很敏感。
两阶段水印（two stages of watermarking，TSW）：将原始图像与水印图像之间的差异可逆地嵌入到水印图像中，以恢复原始图像。通过将比特嵌入图像DCT的AC系数中并将失真嵌入水印图像中。但对几何攻击敏感。对于TSW框架的使用，我们认为可以通过将水印嵌入几何不变域并保证补偿信息足够小来实现。

准备工作

几何攻击

几何变形引起：像素移位，插值导致像素值改变。

泽尼克Zernike矩

In definition, Zernike moments are orthogonal moments based on Zernike polynomials.

图像特征之矩特征本节主要针对Hu矩和Zernike矩作一小结https://mp.weixin.qq.com/s?search_click_id=1160419383886204810-1693894858412-2138666113&__biz=MzU1NDkyMTcxMw==&mid=2247484008&idx=1&sn=339efa27a53763a605b1a87efefb0e1b&chksm=fbdd79f6ccaaf0e0ecbd1cbbad7ab81a130768056d82c6003e4bb1b979a7a37e7a48bec976ba&key=7ec2108f8e296980925bd0ceb0d43643a19ecb38f74d12ca47e3288fcf53ea796d49ed266695b0d79a3844ab3ad982d69099df31f798b6e74e35fc0161e8bde0f156f9d878dcd67c15b801a7790913b0dca4b043abcd9c22c752fe38eb4afda544749e5c4b20c733b907ab296438e605ff895b37c908797b0c76db07d760dfee&ascene=0&uin=MjMwOTg1NzMyOQ%3D%3D&devicetype=Windows+10+x64&version=63090533&lang=zh_CN&countrycode=MV&exportkey=n_ChQIAhIQrgbmeuaEI54TBq3A0hYjABLgAQIE97dBBAEAAAAAABlIJekEgyAAAAAOpnltbLcz9gKNyK89dVj00%2B%2B1e%2BxziTRNQRFSuDmA1ghgqOr%2BNW0YQkAnPOo0eNXi%2FYgnEaxqMFmFIdtKTwyPqkSnzuix214uu9YkRoGWgnZ6jIXO5IZR3C%2FsDW5sS8C1YZk6DvxbhZAVgxTMxFQur%2Fs25poyDwhcMTtOOJDsy2slZa783he7xWYTGjyt8ml4qEWmP0pBceOiNpRN4QqvfrHPsgs1vvi7khFnxgsz%2FC2k0g7tu6Ml9oeOdFjhsa%2B9vairlMw%2FPrFU&acctmode=0&pass_ticket=8MyMl%2FqH%2Ff108tEBivIhVoslL8Hay9aIc4hz9wQk0xI%2BtLaAbrmttFYeoimqqPxW&wx_header=1

距的概念：可以描述全局特征，而且不需要边界限制。

设X为随机变量，c为常数，k为正整数。则量E[(x−c)^k]称为X关于c点的k阶矩。
c=0时，称其为X的k阶原点矩：例如，一阶原点距是期望。
c=E[X]时，称为k阶中心矩。例如，二阶中心距就是方差。

图像中常规距（regular moments）的定义为（上升到二元）：
 

这里距的定义可以理解为，利用图像灰度分布的各阶距来描述图像灰度的分布特性，其中f(x,y)为连续的灰度分布函数（f(x,y)原文中说是图像函数），（p+q）阶距就是以上公式，显然是原点距。
对于数字图像，连续函数就被离散形式所替代：

此时，f(x,y)表示xy处的灰度值，范围为【0,255】

在物理上，m00为灰度质量，m10是x坐标下的质量，m01是y坐标下的质量，那么质心为：

以此类推，那么mpq表示的就是灰度图像相对于原点的分布情况。（太高阶的我也没具体理解）

百度百科-验证（质心的计算）

Hu矩就是定义规格化中心矩，并导出7个不变矩组。

Hu矩在图像描述上有广泛的应用，但是其低阶几何矩与图像整体特征有关，不包含太多的图像细节信息，而高阶几何矩易受噪声影响，因此很难利用几何矩恢复图像。

 regular距的定义是f(x,y)到x^p.y^q上的投影的形式，但是，x^p.y^q不正交，这会导致从距恢复图像困难，所以需要基正交的距，例如Zernike距，此距具有旋转不变性，且便于图像重建，即正交性可以使人们可以分离出每个阶距（信息内容）对重建过程的单独贡献。这些单独贡献简单相加生成重建图像。并且可以通过训练（对比重建图像与原始图像）选择最大阶数，或者加权。

其他矩的概念：k阶矩 (kth-order moment) - 知乎k 阶矩的计算公式（当且仅当下式存在时），其中 F(x) 是随机变量 x 的分布函数： \mathrm{E}\left[X^k\right] =\int_{-\infty}^{+\infty} x^k \mathrm{d}F(x) \\1. 一阶矩 (First order moment)是期望，期望是多次…https://zhuanlan.zhihu.com/p/649317776

引用的Zernike矩原文：Invariant image recognition by Zernike moments | IEEE Journals & Magazine | IEEE Xplore
原文的原文：
Rotation invariant pattern recognition using Zernike moments | IEEE Conference Publication | IEEE Xplore

1979年MichaelReed Teague'所提出Zernike矩：
Captcha

 数学极坐标（r,θ）还有单位元平面中的点Q的笛卡尔坐标(x,y)

Zernike引入一组复多项式，他们在单位圆内部（即x^2+y^2=1）形成一个完全正交集; 因为Zernike多项式具有正交完备性（可以作为基），所以在单位圆中任何图像f(x,y)都可以唯一的被下面的式子展开（这也是正交多项式“正交”的作用，类比xyz坐标，彼此正交，象限中的任何点、函数可以用xyz表示）：

最后θ的定义为p与y轴的夹角，因为原文中定义了：

而且，even为是偶数的意思。那么n与m同奇同偶。

其中，Rnm(p)为径向多项式，定义为（这里没有推导，插个眼！！！）：

Zernike多项式之间具有正交性，其中两个V是任意多项式。最后那个为符号定义，只有a=b时，为1。

Zernike距是图像函数在这些正交基上的投影（拆解）：

其中Anl为Zernike距（投影系数），这也是重建图像的公式。

定义：对于在单位圆外消失的连续图像函数f(x,y)的n阶重复m的距为：

其中V*表示V的共轭,（类似b=ax，那么a=x/b）.
对于数字图像（离散形式）：

且实数展开为：

其中，实部和虚部分别记为Cnl和Snl，则有：

离散化：https://www.cnblogs.com/ronny/p/3985810.html

Zernike距地不变性：

1.旋转不变性对振幅：但是相变了。

2.缩放不变性：但是需要进行归一化。
3.对加类噪声具有鲁棒性：平移、旋转、缩放和噪声操作（有损压缩、低滤波、重采样、重量化、高斯噪声）中保持不变。

前人基于Zernike矩的方法

这些现有的基于Zernike矩的水印方法对内容保持操作提供了很强的鲁棒性，但在可逆性的约束下，由于水印图像与原始图像之间的差异太大而无法可逆嵌入，因此它们都不能用于针对几何变形的无损鲁棒水印。

无损鲁棒水印方法（附代码）

图1展示了三个模块：鲁棒水印的嵌入、补偿信息的可逆嵌入、完整性认证。

1.鲁棒水印的嵌入：首先，从原始图像I中计算泽尼克矩，然后使用私钥随机选择低阶泽尼克矩嵌入水印w，生成中间图像Iw。在嵌入过程中，由水印W引起的失真由三种不同的误差表示。第一个误差来自于量化运算（量化误差dq），第二个误差来自于添加的水印位（水印误差dw），第三种误差来源于对中间图像Iw的舍入操作（舍入误差dr）。

2.补偿信息的可逆嵌入：为了进行完整性认证，对鲁棒水印图像Iw进行舍入运算后Iw1计算哈希值H，并与三个误差放在一起形成补偿信息。最后，将压缩后的补偿信息可逆嵌入到Iwl中，得到最终的水印图像Iw2。

3.完整性认证：在接收端，我们可以提取嵌入的哈希值H来验证图像是否受到攻击，在没有受到攻击的情况下恢复原始图像，或者在受到攻击的情况下直接提取水印进行版权保护。

A.新的水印嵌入策略（水印嵌入到Zernike距中）

图2展示了整个流程：

1.先计算Zernike距

规定m，n

文章中N设定为31（为什么），因为nm同奇同偶，且m最大为31，计算代码如下：

function PQ=zernike_orderlist(ORDER, WITHNEG)

    if nargin<2     %判断函数输入的变量个数，如果小于2，说明WITHNEG没有输入，就置零
        WITHNEG=0;  %如果是1的话，输出的结果包含负数。
    end

    PQ=[];    
    if WITHNEG    
        for p=0:ORDER                         %符合第二条条件
            for q=-p:p                        %此处唯一不同                  
                if mod(abs(p-abs(q)),2)==0    %同奇同偶判断（第一个条件）
                    PQ = [ PQ ; p q ];
                end                    
            end
        end
    else
        for p=0:ORDER
            for q=0:p                  
                if mod(abs(p-q),2)==0
                    PQ = [ PQ ; p q ];
                end                    
            end
        end
    end

计算MxM图像I的单位圆内的Zernike距：

两个步长：

xstep = 2/(img_size(1)-1);  %x方向上步长
ystep = 2/(img_size(2)-1);  %y方向上步长
%因为要在单位圆内计算：
[x,y] = meshgrid(-1:xstep:1,-1:ystep:1); %分别为xy生成-1到1的，且步长为step的矩阵，用于下一步
circle1 = sqrt(x.^2 + y.^2);             %计算距离原点的距离
inside = find(circle1<=1.0001);          %计算位于单位圆内的点
%创建mask矩阵，用来存放该图像中位于单位圆内的点
mask = zeros(img_size);      
mask(inside) = ones(size(inside));  %单位圆内为1，否则为0
Z=zeros(img_size);
Z(inside)=x(inside)+1i*y(inside); %这里没有直接取，而是又另外生成的，与circle1可能边缘有些出入。
p=abs(Z);  %计算模，实际为距离原点的距离
theta=angle(Z);   %最后计算相位 用的是反正切函数

x和y生成的矩阵不同：
图像处理之matlab中meshgrid函数用法详解_meshgrid函数的用法_Hard Coder的博客-CSDN博客1、[X,Y] = meshgrid(x,y) :基于向量x和y中包含的坐标返回二维网格坐标。即A中的每个点对应的是x轴的坐标点，B中的每个点对应的是y轴的坐标点，即有(1,2),(2,2),(3,2),(4,2),(5,2),…,(1,7),(2,7),(3,7),(4,7),(5,7)。3、 [X,Y，Z] = meshgrid(x,y,z),返回由向量x,y,z定义的三维网格坐标，X,Y和Z表示的网格大小为length(x)*length(y)*length(z)。..._meshgrid函数的用法https://blog.csdn.net/qq_44111805/article/details/126273214

Matlab中angle函数内容详解_matlab angle_团子的高速列车的博客-CSDN博客angle函数内容详解作用： 求复数相角，提取相位。用法： p=angle(z);在过去的一周内，由于接触到matlab里的一些函数，需要自己用c代码实现，所以必须去理解函数的作用和操作步骤。这里我用 type +函数名 来获取函数的解释。这里我可以得到函数的具体实现步骤为：而 atan2 函数为反正切函数，atan函数增强版，能确定象限，其函数包含在 ..._matlab anglehttps://blog.csdn.net/cao_jie_xin/article/details/91128887

 接着按照上面提到的公式计算Zernike多项式：对于任一点（x，y）

这里的j是负数 i，原文的定义。

V=zeros(img_size(1),img_size(2),length(PQ));%计算Zernike多项式，有length(PQ)个mn的组合
%按照公式推导，先求R，再求V，（最后求A）
for k=1:length(PQ)
    n=PQ(k,1);
    m=PQ(k,2);
    %一个for循环计算一个Rnm
    a=(n-abs(m))/2;b=(n+abs(m))/2;
    R=zeros(img_size);
    for s=0:a      %这里的小疑问 a与n-|m|/2（公式中的）是否一样？
        num=power(-1,s).*fac(n-s)*(p.^(n-2*s)); %因为p是矩阵，计算出来的是矩阵，对于任意点
        den=fac(s)*fac(b-s)*fac(a-s); %一个实值
        R=R+(num/den); %还是矩阵
    end
    phase=m*theta;
    R(mask==0)=0;
    phase(mask==0)=0;
    V(:,:,k)=(R.*exp(1i*phase));     
end

按照上面计算Zernike距的公式计算：

A=zeros(1,length(PQ)); %求和为一个实值
for k=1:length(PQ)
 n=PQ(k,1);
 A(k)=((n+1)/pi*sum(sum(xstep*ystep*img.*conj(V(:,:,k)))));  %conj计算复数共扼值
end

有有length(PQ)=528个mn的组合，所以计算出528个Zernike距。

2.将比特嵌入选定Zernike距

但是m=4j时，Zernike将偏离正交，所以作为所有符合嵌入条件的Zernike距集合。

in=find(PQ(:,2)>5&PQ(:,1)>2&mod(PQ(:,2),4)~=0);

这里程序中设定为m>5;n>2; m mod 4 != 0; 找到140个。

 假如，需要嵌入的序列为，长度为L，从（私钥）C中随机选取除A00外的距向量并用A00进行归一化，文章用第一个A作为A00。

这里我认为归一化的目的是先归一化到【0,1】：以复数A00为尺，也就是单位，但用T进行放大处理，这样的变换，没有改变其分布，还是可以的。

但是后面的嵌入计算，直接取模，将复数A转变为实数（整数+小数）了。

一般的嵌入和提取：

其中，delta1为量化步长（程序为20），目的是在水印的鲁棒性和不可感知性之间取平衡，Q为量化器：

这是一个向下函数，（如果delta1为1时，则是直接舍弃掉小数）这里delta1为20。

量化器的作用：量化水印(有源代码,绝对可以实现) - 百度文库https://wenku.baidu.com/view/c4a61a236ddb6f1aff00bed5b9f3f90f76c64d94.html?_wkts_=1694750223258&bdQuery=%E6%B0%B4%E5%8D%B0%E5%B5%8C%E5%85%A5%E9%87%8F%E5%8C%96%E5%99%A8%E6%98%AF%E4%BB%80%E4%B9%88

哈哈哈，我的个人白话理解来了：量化操作就是在嵌入时“挖坑”，嵌入就是“埋点土”，但是这个过程我们需要保证水印的保真度（嵌入后的图像与原始图像相近）和鲁棒性（嵌入后可以对抗干扰），因此，量化器就不可以挖的太深，以至于失去保真度，但是也不能太浅，以至于最后提取不出水印。
因为量化误差：是一个实数，需要大量的比特来表示，因此难以嵌入到图像中，所以一般的嵌入不适合。（这里就说挖的太深了）

为了减少量化误差，文章提出了一种新的策略：

嵌入：

其中，（十进制）

下图为水印嵌入的示意图：

a：各个没有嵌入水印之前的计算值。

b:嵌入1之后A所在的位置。

c:嵌入0之后所在的位置。

 提取：

其中， ，，delta为大于0的偶数。

量化误差dq和水印误差dw分别为：

Irw=zeros(img_size);
length_watermark = length(watermark);  %水印长度128
quantified_error = zeros(size(A));    %矩系数绝对值量化损失 dq 528
%水印的嵌入,对zernike矩的幅度值进行量化
for k=1:length_watermark
    A_enc_complex=A(in(k)); %取一个Zernike距
    A_enc_abs = abs(T_multiply*A_enc_complex/A(1)); %归一化 A(1)便是A00
    %  区分出整数部分和小数部分
    [A_enc_int, A_enc_decimal] = Anm_depart(A_enc_abs);
    %选择A的整数部分进行水印嵌入，嵌入后对虚部进行恢复
    q = floor(A_enc_int/step)*step;         %构建量化器q 跟文章不一样 多了一个step
    quantified_error(in(k)) = A_enc_int - q;   %保存整数部分的量化损失
    if watermark(k)==0                  %进行量化水印嵌入
        A_embed=q+step/4;
    else
        A_embed=q+3*step/4;
    end
    Ain=find(abs(A)==abs(A_enc_complex));
    %添加回A小数部分
    A_embed = A_embed + A_enc_decimal;  %跟文章一致，但是文章描述不清晰
    %从加了水印的Zernike距重构为图像形式  这一部分已经到文章第三步了。
    Irw=Irw+(A_embed/A_enc_abs - 1)*(A(in(k))*V(:,:,in(k))+conj(A(in(k)))*V(:,:,Ain(1)));   
end

%将高频成分和带有水印的低频成分加入，得到水印图像Iw
Iw=img+Irw;
%对Iw进行round处理，保证Iw为整数
Iw1=round(Iw);

3.可逆Zernike变换

原始的Zernike距为，归一化之后的距为，
利用振幅生成水印版的Zernike距。
此为加了水印后的Zernike矩向量。

重构水印：

划线的就是水印

原理：

 加上水印：

加完可能会超过255，需要标准化舍入一下。

watermarked_error=watermark;
%溢出损失，程序实现与文章有点区别
overflow_error = double(Iw1) - double(uint8(Iw1)); %计算上，没有溢出损失

4.舍入误差的计算（第三个误差）

Usually, the distortion is represented as the difference between the inter mediate watermarked image I w and the robust watermarked image I w 1 .

与上一个过程相同，计算Zernike距，但是带入的是加水印的图像：

恢复原始的Zernike矩

水印影响，这里看清是A-Aw，所以是水印的负值。

恢复原始图像

但还是存在差异：

%----------------3.处理取整失真----------------
%zernike正变换，得到Aw', V
[Aw]=zernike_A(Iw1, V, PQ);
Ir=zeros(img_size);
% 选择合适的zernike矩
for k=1:length_watermark
    Aw_enc_complex=Aw(in(k)); %之前的选择
    Aw_enc_chosen = abs(T_multiply*Aw_enc_complex/Aw(1)); %归一化
     %选择A的整数进行水印提取，提取后对小数部分和量化损失进行恢复
    [Aw_enc_int, Aw_dec_decimal] = Anm_depart_decode(Aw_enc_chosen,step);
    if watermark(k)==1
        A_enc_recover = (Aw_enc_int - 3/4*step + quantified_error(in(k))) + Aw_dec_decimal;
    else
        A_enc_recover = (Aw_enc_int - 1/4*step + quantified_error(in(k))) + Aw_dec_decimal;
    end  
    Aw_enc_in=find(abs(Aw)==abs(Aw_enc_complex));
    Ir = Ir+(A_enc_recover/Aw_enc_chosen -1 )*(Aw(in(k))*V(:,:,in(k)) + conj(Aw(in(k)))*V(:,:,Aw_enc_in(1)));  
end
% 加上高频成分重构出原始图像
img_pre_recovered = Iw1 + Ir;
img_pre_recovered = round(img_pre_recovered);
%计算取整损失 rounded_loss
rounded_error = double(img) - double(img_pre_recovered);
%得到0-255的水印图像
Iw1_uint8 = uint8(Iw1);
%生成Iw1_uint8的哈希值 主要用于下下步验证
hash_value=hash(char(Iw1_uint8),'sha256'); 

%存在的失真有：量化损失 dq，取整损失 dr, 溢出损失 do 水印损失dw

5.补偿信息的可逆嵌入

噪声影响后，补偿信息：三个误差以及Iw1的hash值。以可逆水印的方式嵌入（这里埋个坑，因为代码太长了，而且是基于提出的文献，这里就不展开了）

B.完整性认证

主要是比较前后两个图像的hash值。

C.提取水印和恢复图像（前提是不受攻击）

1）补偿信息的提取 2）水印的提取 3）恢复图像

D.在攻击在的水印提取

整个过程已经了解，通了一遍，还差结果部分，在此就不写了。