0引言随着互联网的迅速发展，人们与同事、朋友、家人等分享信息变得更加方便，可以便捷地将图像、视频等信息通过网络传播。然而，在享受这种便捷的同时，图像、视频等载体中的个人信息以及分享者的合法权益都会受到威胁。为解决这一问题，数字媒体的安全与保护成为重要的研究课题。通常，数字媒体包括音频，图像和视频；其中，图像是传递信息的主要途径。通过将图像看做2维信号，传统数字信号处理领域的许多技术都可以用来处理图像。作为离散傅里叶变换(discrete fourier transform，DFT)的推广，离散分数阶傅里叶变换(discrete fractional fourier transform，DFrFT) 广泛应用于信号和图像处理领域，如信号分类(詹宏锋，2020)、图像加密(Hsue，2019；辛怡等，2008)、图像融合(Saxena和Sharma，2017)和数字水印(Tsai和Hsue，2015)。DFT的其他应用可参阅Tao等人(2008)。事实上，通常图像的像素是0~255之间的实整数。然而，DFrFT变换是一种复数变换，图像的DFrFT变换将实值图像域转化为复数域。这会导致变换系数的存储量和计算复杂度的增加。在某些情况下，将实值数据转换为变换域的实值数据系数是非常有用的。因此，诸多学者研究了实离散分数余弦和正弦变换(Venturini和Duhamel，2004)。然而，在Venturini和Duhamel(2004)的工作中，实离散分数余弦和正弦变换矩阵被限制为长度N=4$m$($m$=1, 2, 3…)。Hsue和Chang(2015)给出了构造任意长度实离散分数傅里叶变换的一般方法，并定义了实离散分数傅里叶变换(real discrete fractional Fourier transform，RDFrFT)。Yap等人(2003)提出了Krawtchouk变换(Krawtchouk transform，KT)，该变换在图像分析领域具有诸多优点，得到广泛研究。根据KT变换矩阵的对称性，Abdulhussain等人(2018)研究了KT的快速计算方法。KT变换的理论研究保证了其在图像识别和图像表示中的广泛应用(Xiao等，2016；Abdulhussain等，2019)。作为传统KT的推广，Liu等人(2017)提出了离散分数阶Krawtchouk变换(discrete fractional KT, DFrKT)，并将DFrKT应用于图像水印。结果显示，DFrKT在图像水印应用中具有较好的表现。本文将Krawtchouk变换推广到实离散分数Krawtchouk变换(real DFrKT, RDFrKT)，并给出了RDFrKT的一些性质。为了验证RDFrKT的实用性，RDFrKT被应用于数字图像水印算法的设计。最后，实验结果表明，基于RDFrKT的水印方案具有较好的水印鲁棒性。1相关理论一维$N \times N$大小的Krawtchouk变换矩阵$\boldsymbol{K}$(Liu等，2017)定义为 1 $\boldsymbol{K}_{n, x}=K_{n}(x ; p, N-1), 0 \leqslant n, x \leqslant N-1$ 式中，${K_n}\left({x;p, N - 1} \right)$是$n$阶加权Krawtchouk多项式，定义为 2 $K_{n}(x ; p, N-1)=r_{n}(x ; p, N-1) \sqrt{\frac{w(x ; p, N-1)}{\rho(n ; p, N-1)}}$ 并且 3 $w(x ; p, N-1)=\left(\begin{array}{c}N-1 \\x\end{array}\right) p^{x}(1-p)^{N-1-x}$ 4 $\rho(n ; p, N-1)=\left(\frac{p-1}{p}\right)^{n} \frac{n !}{(-N+1)_{n}}$ 式中，$\left({\begin{array}{*{20}{c}} {N - 1}\\ x \end{array}} \right)$为二项分布系数，${r_n}\left({x;p, N - 1} \right)$是Krawtchouk多项式，定义为 5 $\begin{gathered}r_{n}(x ; p, N-1)= \\{ }_{2} F_{1}\left(-n,-x ;-N+1 ; \frac{1}{p}\right), p \in(0,1)\end{gathered}$ 式中，${}_2F$表示超几何函数。矩阵$\boldsymbol{K}$的特征值分解(Xiao等，2016；肖斌等，2018)表示为 6 $\boldsymbol{K}=\boldsymbol{V} \boldsymbol{D} \boldsymbol{V}^{\mathrm{T}}$ 式中，矩阵$\boldsymbol{V}$表示$\boldsymbol{K}$的正交特征向量，对角矩阵$\boldsymbol{D}$表示特征值矩阵，其对角元素为$\boldsymbol{K}$的特征值。如果$\boldsymbol{N}$为偶数，$\boldsymbol{D}$的特征值排列为 7 $\boldsymbol{D}=\left[\begin{array}{lllll}1 & & & & \\& -1 & & & \\& & 1 & & \\& & & \ddots & \\& & & & -1\end{array}\right]$ 如果$\boldsymbol{N}$为奇数，则$\boldsymbol{D}$的特征值排列为 8 $\boldsymbol{D}=\left[\begin{array}{lllll}1 & & & & \\& -1 & & & \\& & 1 & & \\& & & \ddots & \\& & & & 1\end{array}\right]$ 此外，矩阵$\boldsymbol{K}$的特征值分解还可以表示为 9 $\boldsymbol{K}=\sum\limits_{k=0}^{N-1} d_{k} \boldsymbol{v}_{k} \boldsymbol{v}_{k}^{\mathrm{T}}$ 式中，${d_0}, {d_1}, \cdots, {d_{N - 1}}$是$\boldsymbol{K}$的特征值，其值为$\boldsymbol{D}$的主对角元素。${\boldsymbol{v}_0}, {\boldsymbol{v}_1}, \cdots, {\boldsymbol{v}_{N - 1}}$是分别对应特征值${d_k}\left({k = 0, 1, \cdots, N - 1} \right)$的$\boldsymbol{K}$的一组正交特征向量。根据式(9)，通过对$\boldsymbol{K}$的特征值取分数阶幂，可以构造离散分数阶Krawtchouk变换(DFrKT)矩阵，$\alpha $阶DFrKT变换矩阵${\boldsymbol{K}^\alpha }$定义为 10 $\boldsymbol{K}^{\alpha}=\sum\limits_{k=0}^{N-1} \mathrm{e}^{-\mathrm{j} k {\rm{ \mathsf{ π} }} \alpha} \boldsymbol{v}_{k} \boldsymbol{v}_{k}^{\mathrm{T}}$ 式中，分数阶$\alpha $取值范围为[0, 1]，${\boldsymbol{v}_0}, {\boldsymbol{v}_1}, \cdots, {\boldsymbol{v}_{N - 1}}$是$K$的一组正交特征向量。那么图像$\boldsymbol{f}$的($\alpha $, $\beta $) 阶2维DFrKT可以表示为 11 $\boldsymbol{T}^{\alpha, \beta}=\boldsymbol{K}^{\alpha} \boldsymbol{f} \boldsymbol{K}^{\beta}$ 从式(11)可以发现，DFrKT的变换矩阵${\boldsymbol{K}^\alpha }$和${\boldsymbol{K}^\beta }$都是复值矩阵，因此，图像的DFrKT变换将实值图像域变换到了复数域。相应的DFrKT逆变换可表示为 12 $\boldsymbol{f}=\boldsymbol{K}^{-\alpha} \boldsymbol{T}^{\alpha, \beta} \boldsymbol{K}^{-\beta}$ 从式(12)可以发现，如果变换系数${\boldsymbol{T}^{\alpha, \beta }}$被误差干扰，重构图像$\boldsymbol{f}$将是复值函数。如果使用被误差干扰的扰动系数来重建图像，得到的图像域是个复数域，需要额外的存储空间来储存重建图像。2实离散分数阶Krawtchouk变换2.1RDFrKT的定义根据DFrKT的定义可知，DFrKT变换矩阵的构造是通过对Krawtchouk变换矩阵$\boldsymbol{K}$进行特征值分解后，对特征值取分数次幂得到的。为了构造变换域为实数的DFrFT变换来处理实值信号，Hsue和Chang(2015)给出了一种根据变换矩阵的特征值分解来构造实离散分数阶Fourier变换的方法。本文将该方法推广到RDFrKT的定义。定义大小为$N \times N$的$\alpha $阶实离散分数Krawtchouk变换矩阵$\boldsymbol{K}_r^\alpha $为 13 $\boldsymbol{K}_{r}^{\alpha}=\overline{\boldsymbol{V}} \overline{\boldsymbol{D}}{}^{\alpha} \overline{\boldsymbol{V}}{}^{\mathrm{T}}$ 式中，$\boldsymbol{\bar V}$的每列是矩阵$\boldsymbol{K}$的一组正交特征向量，其第$m$列为${\boldsymbol{\bar v}_m}$。${\boldsymbol{\bar D}^\alpha }$是分块对角矩阵，表示为 14 $\overline{\boldsymbol{D}}{}^{\alpha}=\left[\begin{array}{lllll}\overline{\boldsymbol{D}}{}_{0,0}^{\alpha} & & & & \\& \overline{\boldsymbol{D}}{}_{1,1}^{\alpha} & & & \\& & \overline{\boldsymbol{D}}{}_{2,2}^{\alpha} & & \\& & & \ddots & \\& & & & \overline{\boldsymbol{D}}{}_{M-1, M-1}^{\alpha}\end{array}\right]$ 式中，$\boldsymbol{\bar D}_{m, m}^\alpha $($m$=0, 2, …, $M$-1)的构造与${\boldsymbol{\bar D}^\alpha }$矩阵的大小$M$有关，即1) 当$M = 4 \times \left\lfloor {N/4} \right\rfloor $，$\boldsymbol{\bar D}_{m, m}^\alpha $ ($m$=0, 2, …, $M$-1)是Givens矩阵，定义为 15 $\overline{\boldsymbol{D}}{}_{m, m}^{\alpha}=\left[\begin{array}{cc}\cos (\alpha m {\rm{ \mathsf{ π} }}) & \sin (\alpha m {\rm{ \mathsf{ π} }}) \\-\sin (\alpha m {\rm{ \mathsf{ π} }}) & \cos (\alpha m {\rm{ \mathsf{ π} }})\end{array}\right]$ 2) 当$M = 4 \times \left\lfloor {N/4} \right\rfloor $+1，$\boldsymbol{\bar D}_{m, m}^\alpha $($m$=0, 2, …, $M$-1)定义为式(15)，且$\boldsymbol{\bar D}_{M, M}^\alpha = \boldsymbol{I}$, $\boldsymbol{I}$为单位矩阵。3) 当$M = 4 \times \left\lfloor {N/4} \right\rfloor $+2，$\boldsymbol{\bar D}_{m, m}^\alpha $($m$=0, 2, …, $M$-2)定义为式(15)，且$\boldsymbol{\bar D}_{M - 1, M - 1}^\alpha = \boldsymbol{I}$，$\bar D_{M, M}^\alpha = - \boldsymbol{I}$。4) 当$M = 4 \times \left\lfloor {N/4} \right\rfloor $+3，$\boldsymbol{\bar D}_{m, m}^\alpha $($m$=0, 2, …, $M$-1)定义为式(15)，且$\boldsymbol{\bar D}_{M, M}^\alpha = - \boldsymbol{I}$。然后构造特征向量vm，对于N4，vm($m$=0, 1, 2, …, $N$-1), 构造规则如下：1) 若$m$=4k或4k+1, (k=0, 1, …, $\left\lfloor {N/4} \right\rfloor - 1$), vm是$\boldsymbol{K}$对应特征值1的正交向量；若$m$=4k+2或4k+3 (k=0, 1, …, $\left\lfloor {N/4} \right\rfloor - 1$)，vm是$\boldsymbol{K}$对应特征值-1的正交向量。2) 若$N = 4 \times \left\lfloor {N/4} \right\rfloor$ +1，则${\boldsymbol{\bar v}_N}$为$\boldsymbol{K}$对应特征值1的正交向量。3) 若$N = 4 \times \left\lfloor {N/4} \right\rfloor$+2，则${\boldsymbol{\bar v}_{N - 2}}$为$\boldsymbol{K}$对应特征值1的正交向量，${\boldsymbol{\bar v}_{N - 1}}$为$\boldsymbol{K}$对应特征值-1的正交向量。4) 若$N = 4 \times \left\lfloor {N/4} \right\rfloor$+3，则${\boldsymbol{\bar v}_{N - 3}}$, ${\boldsymbol{\bar v}_{N - 2}}$为$\boldsymbol{K}$对应特征值1的正交向量，${\boldsymbol{\bar v}_{N - 1}}$为$\boldsymbol{K}$对应特征值-1的正交向量。5) 若$N = 4 \times \left\lfloor {N/4} \right\rfloor$+4，则${\boldsymbol{\bar v}_{N - 3}}$, ${\boldsymbol{\bar v}_{N - 4}}$为$\boldsymbol{K}$对应特征值1的正交向量，${\boldsymbol{\bar v}_{N - 2}}$，${\boldsymbol{\bar v}_{N - 1}}$为$\boldsymbol{K}$对应特征值-1的正交向量。值得注意的是，当N≤4时，可以使用规则2)—4)生成${\boldsymbol{\bar v}_m}$。根据实离散分数Krawtchouk变换矩阵的定义式(13)，可以定义作用于图像$\boldsymbol{f}$的2维实离散分数Krawtchouk变换(RDFrKT) $\boldsymbol{T}_r^{\alpha, \beta }$为 16 $\boldsymbol{T}_{r}^{\alpha, \beta}=\boldsymbol{K}_{r}^{\alpha} \boldsymbol{f} \boldsymbol{K}_{r}^{\beta}$ 式中，$\boldsymbol{K}_r^\alpha $和$\boldsymbol{K}_r^\beta $分别是分数阶为$\alpha, \beta $的RDFrKT变换矩阵。根据DFrKT变换矩阵(13)的定义，可以看出，$\boldsymbol{K}_r^\alpha $的逆可以通过取分数阶数$\alpha $的负数得到，即 17 $\boldsymbol{K}_{r}^{\alpha} \boldsymbol{K}_{r}^{-\alpha}=\boldsymbol{I}$ 根据式(17)，2维逆RDFrKT可以定义为 18 $\boldsymbol{f}=\boldsymbol{K}_{r}^{-\alpha} \boldsymbol{T}_{r}^{\alpha, \beta} \boldsymbol{K}_{r}^{-\beta}$ 2.2RDFrKT变换矩阵的性质通常而言，构造的分数阶变换矩阵应满足3个条件：一致性、阶数可加性和幺正性。对所提出的RDFrKT而言，传统的KT变换是RDFrKT的一个特例，一致性条件满足；并且它具有阶数可加性和幺正性。RDFrKT的性质解释如下：性质1：如果RDFrKT变换矩阵$\boldsymbol{K}_r^\alpha $的分数阶$\alpha $=0，则$\boldsymbol{K}_r^0 = \boldsymbol{I}$。性质2：一致性。如果RDFrKT变换矩阵$\boldsymbol{K}_r^\alpha $的分数阶$\alpha $=1，则$\boldsymbol{K}_r^1 = \boldsymbol{K}$。这表明当分数阶为1时，RDFrKT将变为Krawtchouk变换。性质3：阶数可加性，即 19 $\boldsymbol{K}_{r}^{\alpha} \boldsymbol{K}_{r}^{\beta}=\boldsymbol{K}_{r}^{\alpha+\beta}$ 证明：由式(13)可得 20 $\boldsymbol{K}_{r}^{\alpha} \boldsymbol{K}_{r}^{\beta}=\overline{\boldsymbol{V}} \overline{\boldsymbol{D}}{}^{\alpha} \overline{\boldsymbol{V}}{}^{\mathrm{T}} \overline{\boldsymbol{V}} \overline{\boldsymbol{D}}{}^{\beta} \overline{\boldsymbol{V}}{}^{\mathrm{T}}$ 然后，利用$\boldsymbol{\bar V}$正交性${\boldsymbol{\bar V}^{\rm{T}}}\boldsymbol{\bar V} = \boldsymbol{I}$，式(20)可以进一步表示为 21 $\boldsymbol{K}_{r}^{\alpha} \boldsymbol{K}_{r}^{\beta}=\overline{\boldsymbol{V}} \overline{\boldsymbol{D}}{}^{\alpha} \overline{\boldsymbol{D}}{}^{\beta} \overline{\boldsymbol{V}}{}^{\mathrm{T}}=\overline{\boldsymbol{V}} \overline{\boldsymbol{D}}{}^{\alpha+\beta} \overline{\boldsymbol{V}}{}^{\mathrm{T}}=\boldsymbol{K}_{r}^{\alpha+\beta}$ 即证。性质4：幺正性，即 22 $\boldsymbol{K}_{r}^{-\alpha}=\left(\boldsymbol{K}_{r}^{\alpha}\right)^{-1}$ 性质4的证明可通过使$\beta $=-$\alpha $和使用性质1来实现。3在图像水印中的应用相对于复变换域，信号的实变换域具有在计算机中存储时占用空间少、在数字通信信道占用带宽少的优点，实值变换已用于图像加密(Kang等，2019；Lang，2012；Zhou等，2012)和图像分类(Zhang等，2018)等图像处理应用中。作为RDFrKT的一个应用，基于所提出的RDFrKT，本文提出一种结合RDFrKT和奇异值分解(singular value decomposition, SVD)数字图像水印方案。基于SVD的水印方案一般通过两种方式嵌入水印：通过修改奇异值嵌入水印(Mohammad等，2008；Run等，2012；Pei和Liu，2008；Guo和Prasetyo，2014)和将水印嵌入奇异向量(Jia，2014；Hu等，2020)。而且，为了增强水印的鲁棒性，可以将水印嵌入到图像DCT(discrete cosine transform)变换域的SVD中(Tian等，2020；Singh和Singh，2017)。实值变换将实值的图像域映射到实值的变换域，对实值的变换域处理后逆变换域仍然是实值。这样，对图像的实变换域修改后逆变换得到的仍然是实的图像域。因此，在图像的实变换域嵌入水印信息是理想的选择。本文通过图像RDFrKT的变换域嵌入二值水印信息，将水印信息嵌入到RDFrKT变换域的奇异值分解后的奇异向量中。水印方案分为两个阶段：水印嵌入阶段和水印提取阶段。图 1和图 2分别直观地给出了水印的嵌入和提取流程。 图1 水印嵌入流程 The proposed watermark embedding procedureFig 1 图2 水印提取流程 The proposed watermark extraction procedureFig 23.1水印嵌入流程考虑$N \times N$大小的原始灰度图像$\boldsymbol{A}$和$m$×$m$大小的水印二值图像$\boldsymbol{W}$。水印的嵌入流程包括以下几个步骤：1) 将原始图像$\boldsymbol{A}$分成4×4大小的不重叠的图像块。2) 为了加强水印的安全性，随机选取m2(m≤N/4) 个图像块。图像块的选取位置可以看做秘钥$\boldsymbol{KE}{\boldsymbol{Y}_1}$。3) 对每个块执行RDFrKT，得到变换域系数$\boldsymbol{X}$。4) 对每个块的RDFrKT系数执行SVD分解，即 23 $\boldsymbol{X}=\boldsymbol{U} \boldsymbol{S} \boldsymbol{V}^{\mathrm{T}}$ 式中, $\boldsymbol{U}$为$\boldsymbol{X}$的左奇异向量矩阵，$\boldsymbol{V}$为右奇异向量矩阵，$\boldsymbol{S}$为奇异值矩阵。5) 通过修改左奇异向量矩阵第1个奇异向量的第2个元素u21和第3个元素u31的差值来嵌入1 bit水印w。设$g = {u_{21}} - {u_{31}}$，水印块使用以下扰动调制的方法(Hu等，2020)通过修改$g$来嵌入水印，即 24 $g^{*}= \begin{cases}\gamma_{21} & g\rho \text { 且 }\left|\gamma_{21}-g\right|\left|\gamma_{1}-g\right| \\ \gamma_{22} & g-\rho \text { 且 }\left|\gamma_{22}-g\right|\left|\gamma_{1}-g\right| \\ \gamma_{1} & \text {其他 }\end{cases}$ 式中，$g^*$是对$g$嵌入水印后的值，${\gamma _{21}}$和${\gamma _{22}}$表示为 25 $\gamma_{1}= \begin{cases}\min (\psi, \max (g, \delta)) & w=1 \\ \max (-\psi, \min (g,-\delta)) & w=0\end{cases}$ 26 $\begin{gathered}\gamma_{21}= \\\begin{cases}\max \left(0,\left\lfloor {\frac{g-\eta}{2 D}} \right\rfloor \times 2 D+D+\eta\right) & w=1 \text { 且 } g\rho \\{\left\lfloor {\frac{g-\eta}{2 D}+0.5} \right\rfloor \times 2 D+\eta} & w=0 \text { 且 } g\rho\end{cases}\end{gathered}$ 27 $\begin{gathered}\gamma_{22}= \\\begin{cases}{\left\lceil {\frac{g+\eta}{2 D}-0.5} \right\rceil \times 2 D-\eta} & w=1 \text { 且 } g-\rho \\\min \left(0,\left\lceil {\frac{g+\eta}{2 D}} \right\rceil \times 2 D-D-\eta\right) & w=0 \text { 且 } g-\rho\end{cases}\end{gathered}$ 式中，max()和min()分别表示函数的最大值和最小值。式(24)—(27)给出的参数$\rho, \eta, D$，$\delta = 0.5D, \psi = D, \rho = 1.5D, \eta = 2D$，调整量化步长$D$来平衡水印的鲁棒性和不可见性。6) 在修改$g$嵌入水印之后，可以得到修改后的奇异向量元素${u_{21}}$和${u_{31}}$为 28 $\left\{\begin{array}{l}u_{21}^{*}=u_{21} \\u_{31}^{*}=u_{21}-g^{*}\end{array}\right.$ 7) 得到修改后的RDFrKT系数$\boldsymbol{X}^*$为 29 $\boldsymbol{X}^{*}=\boldsymbol{U}^{*} \boldsymbol{S} \boldsymbol{V}^{\mathrm{T}}$ 式中，$\boldsymbol{U}^*$为嵌入水印后图像块的左奇异向量矩阵。8) 对每个块的修改后的RDFrKT系数$\boldsymbol{X}^*$执行逆RDFrKT来获得加水印后的图像$\boldsymbol{A}^*$。3.2水印提取流程假如接收到一幅图像$\boldsymbol{A'}$，必须拥有秘钥$\boldsymbol{KE}{\boldsymbol{Y}_1}$，才能正确地提取出水印信息。详细的水印提取步骤如下：1) 将接收到的图像$\boldsymbol{A'}$分成4×4大小的不重叠的图像块。2) 使用$\boldsymbol{KE}{\boldsymbol{Y}_1}$定位正确水印的图像块。3) 对每个块Block(i, j)进行RDFrKT来获得变换域系数$\boldsymbol{X'}$，接着对$\boldsymbol{X'}$进行SVD分解，得到 30 $\boldsymbol{X}^{\prime}=\boldsymbol{U}^{\prime} \boldsymbol{S}^{\prime} \boldsymbol{V}^{\prime \mathrm{T}}$ 式中，$\boldsymbol{U'}$为$\boldsymbol{X'}$的左奇异向量矩阵，$\boldsymbol{V'}$为右奇异向量矩阵，$\boldsymbol{S'}$为特征值矩阵。4) 从左奇异向量矩阵的第2和第3个元素${u'_{21}}$和${u'_{31}}$之间的差$g'$中提取水印比特，即 31 $\begin{gathered}w^{\prime}= \\\begin{cases} 0 & \left|g^{\prime}\right|\rho \text { 且 } g^{\prime} \geqslant 0 \\1 & \left|g^{\prime}\right|\rho \text { 且 } g^{\prime}0 \\ \bmod \left(\left\lfloor{\frac{g^{\prime}-\eta}{D}+0.5} \right\rfloor, 2\right) & g^{\prime}\rho \\ \bmod \left(\left\lceil {\frac{g^{\prime}+\eta}{D}-0.5}\right\rceil-1,2\right) & g^{\prime}-\rho\end{cases}\end{gathered}$ 式中，$\bmod \left({ \cdot, 2} \right)$表示取模运算，$w'$为一个图像块提取的1 bit水印。最后，组合所有图像块提取的水印信息可以用来确定图像的所有权或验证版权。4实验结果及分析为了验证提出的基于RDFrKT的水印算法的有效性，实验采用了Granada大学计算机视觉组的图像数据库(http://decsai.ugr.es/cvg/dbimagenes/)中的96幅大小为256×256像素的灰度图像为测试图像。采用MPEG-7二值图像形状数据库(http://imageprocessingplace.com/root_files_V3/image_databases.htm)中的64×64像素“Deer”二值图像为水印。图 3给出了几种测试图像和水印图像实例。实验中，水印图像嵌入到测试图像后生成嵌入水印的载体图像。然后，载体图像经过滤波、噪声等图像攻击后，测试从受攻击的载体图像提取水印信息的有效性，从而验证所提算法的有效性。为了验证所提水印算法的鲁棒性能，将本文方法与基于DCT、KT的方法和基于Hu的SVD方法(Hu等，2020)进行了比较。同时，为了保证比较的公平性，本文水印算法同时考虑了修改奇异向量时产生的移位、符号矫正、正交恢复等预处理操作(Hu等，2020)。 图3 测试图像数据库和二值水印图像的若干样本 Some samples of the test image database and binary watermark imageFig 3为了度量水印在载体图像中的不可见性，实验采用了峰值信噪比(peak signal to noise ratio，PSNR)来衡量，其定义为(Run等，2012) 32 $P S N R=10 \lg \frac{255^{2}}{M S E}$ 式中，$MSE$表示水印图像${\boldsymbol{A}^*}$和原始图像${\boldsymbol{A}^*}$之间的均方误差，定义为 33 $M S E=\frac{1}{N N} \sum\limits_{i=1}^{N} \sum\limits_{j=1}^{N}\left(\boldsymbol{A}^{*}-\boldsymbol{A}\right)$ 为了保证嵌入水印后载体图像的视觉质量，实验对水印参数进行了调整，使测试图像的平均PSNR为37.58 dB。然后选择量化步长$D$=0.025。根据实验测试，为了得到鲁棒的水印算法，RDFrKT变换中分数阶阶数的选取为$\alpha $=0.6，$\beta $=0.4。水印鲁棒性用误码率(bit error rate，BER)来衡量，定义为(Tian等，2020) 34 $B E R=\frac{\sum\limits_{i=1}^{l} \sum\limits_{j=1}^{l}\left|\boldsymbol{W}^{*}-\boldsymbol{W}\right|}{l \times l}$ 式中，${\boldsymbol{W}^*}$是提取的水印图像，$\boldsymbol{W}$是尺寸为$l \times l$的原始二值水印。首先以Lena图像作为原始图像，将“Deer”水印嵌入到原始图像中。然后，嵌入水印的载体图像受到各种类型的水印攻击，如中值滤波、均值滤波、高斯滤波、JPEG压缩、图像缩放、椒盐噪声和高斯噪声。然后从受攻击后的载体图像中提取水印。图 4给出了一些提取的水印及其在不同攻击情况下对应的BER值。可以看出，在大多数攻击情况下，基于RDFrKT方法提取的水印比其他方法具有更低的BER值。这表明基于RDFrKT的方法具有更好的性能。 图4 不同攻击的提取水印及其误码率(BER) Some example of extracted watermarks and their BER values of different attacksFig 4 ((a)ours; (b)DCT; (c)KT; (d)Hu et al.(2020)) 图像库中的96幅图像为测试图像，将水印分别嵌入其中得到嵌入水印的载体图像。然后，所有的载体图像分别在中值滤波、均值滤波、高斯滤波、JPEG压缩、图像缩放、椒盐噪声和高斯噪声等各种水印攻击后提取嵌入信息。统计了所有测试图像的平均BER值来验证本文所提算法的鲁棒性能。图 5给出了各种类型攻击下，平均BER值随攻击强度的变化情况，并将本文算法与基于DCT的水印算法、基于KT的水印算法、Hu等人(2020) SVD水印算法进行了比较。结果表明，基于DFrKT的水印算法在大多数情况下都具有更好的鲁棒性，特别是对中值滤波、均值滤波、高斯滤波和JPEG压缩这几类图像攻击。 图5 各种攻击类型下平均BER值随攻击强度的变换情况 Average BER value variation with attack parameter changes of different attack situationsFig 5 ((a) median filter; (b) average filter; (c) Gaussian filter; (d) JPEG compression; (e) image recalling; (f) salt and peppers noise; (g) Gaussian noise) 5结论本文提出了实离散分数Krawtchouk变换(RDFrKT)。与现有的离散分数Krawtchouk变换(DFrKT)相比，RDFrKT可以将实数值的灰度图像域变换到实数值的变换域。于是，RDFrKT变换域具有和图像相同的存储空间；并且，图像的RDFrKT变换域进行处理后进行反变换，得到的图像空间仍然是实数域。因此，RDFrKT不仅保持了传统Krawtchouk变换是实变换的特点，并且具有可随分数阶阶数调整的变换域，在图像处理中具有广泛的应用前景。作为RDFrKT的一个应用，该变换已用于数字图像水印算法的设计。结果表明，基于RDFrKT变换的水印方案，对图像中值滤波、均值滤波、高斯滤波、JPEG压缩和缩放等多种信号处理攻击具有较好的鲁棒性。本文设计的水印算法结合了RDFrKT和SVD，在图像的奇异向量中嵌入了水印信息。在未来的工作中，RDFrKT可以与SVD结合，将水印信息嵌入到奇异值来探索水印算法的性能。此外，所提出的变换有望应用于图像加密和其他水印应用中。一种思路是可以将辛怡等人(2008)的实分数阶Fourier变换域图像加密算法扩展到RDFrKT变换域。