信息隐藏与数字水印实验3-DCT 域图像水印

一、 实验目的
了解频域水印的特点，掌握基于DCT系数关系的图像水印算法原理，设计并实现一种基于DCT域的图像水印算法，嵌入二值图像水印信息，掌握水印图像的归一化函数的计算方法，并对携秘图像进行攻击，提取攻击后的水印二值图像，计算NC的值。

二、 实验环境
(1) Windows或Linux换作系统
(2) Python3 环境
(3) Python的 opencv-python、 numpy、 matplotlib、 pycrypto 库
(4)原始载体图像文件和二值水印图像文件

三、 实验要求
（1）利用文件夹中的彩色图片lena.bmp图片在B颜色通道嵌入秘密信息，这个秘密信息为10*10的二值图像水印bupt10.bmp。嵌入后的图像为lenastegoB.bmp图像。在一行两列中输出lena.bmp和lenastegoB.bmp图像。输出图像的时候需要在标题处增加本人学号。
（2）计算嵌入秘密信息后lenastegoB.bmp的峰值信噪比。输出峰值信噪比，要求峰值信噪比精确到小数点后三位。
（3）在没有任何攻击的情况下从lenastegoB.bmp中图像提取水印二值图像bupt1.bmp，计算bupt10.bmp和bupt1.bmp的归一化函数NC，在没有任何攻击的情况下，提取出来的水印信息和原始的水印信息完全相同，归一化函数的值为1。归一化的计算方法参考Wav音频的LSB信息隐藏。
（4）在lenastegoB.bmp图像中添加白噪声，将图像另存为lenastegoBl.bmp。从
添加白噪声后的携密图像lenastegoB1.bmp中提取水印二值图像bupt2.bmp，计算bupt10.bmp和bup2.bmp的归一化函数NC，因为携密图像遭受了高斯噪声攻击，提取出来的水印信息和原始的水印信息完全相同，归一化函数的值肯定不为1。计算这个NC值并输出。要求NC的值精确到小数点后三位。

四、 实验步骤和结果

作业：

（1）利用文件夹中的彩色图片lena.bmp图片在B颜色通道嵌入秘密信息，这个秘密信息为10*10的二值图像水印bupt10.bmp。嵌入后的图像为lenastegoB.bmp图像。在一行两列中输出lena.bmp和lenastegoB.bmp图像。输出图像的时候需要在标题处增加本人学号。

实验步骤：

（1）以灰度图模式读取水印图像，将其从二维矩阵转为一维并二值化，随后读取载体图像lena.bmp。
（2）随后取图片的蓝色分量隐藏，对图像分块并进行信息嵌入，将分块合并为完整的图像蓝色分量，并与图像的绿色，红色分量合并。
（3）最后计算峰值信噪比，保存嵌入水印的载体图像，并展示原图像和嵌入图像。

block.py代码

import numpy as np


def img_to_blocks(img, block_shape):
    height, width = img.shape[:2]
    block_height, block_width = block_shape
    shape = (height // block_height, width // block_width, block_height, block_width)
    strides = img.itemsize * np.array([width * block_height, block_width, width, 1])
    img_blocks = np.lib.stride_tricks.as_strided(img, shape, strides).astype('float32')
    img_blocks = np.reshape(img_blocks, (shape[0] * shape[1], block_height, block_width))
    return img_blocks


def blocks_to_img(img_blocks, img_shape):
    height, width = img_shape[:2]
    block_height, block_width = img_blocks.shape[-2:]
    shape = (height // block_height, width // block_width, block_height, block_width)
    img_blocks = np.reshape(img_blocks, shape)

    lines = []
    for line in img_blocks:
        lines.append(np.concatenate(line, axis=1))
    img = np.concatenate(lines, axis=0)
    return img

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25

代码：

# 导入图像处理库
import cv2
# 导入数学计算库
import numpy as np
# 导入绘图库
import matplotlib.pyplot as plt
# 从block.py导入图像分块的函数、分块合并的函数
from block import img_to_blocks, blocks_to_img


# 计算峰值信噪比的函数
def PSNR(template, img):
    mse = np.mean((template / 255. - img / 255.)**2)
    if mse < 1.0e-10:
        return 100
    PIXEL_MAX = 1
    return 20 * np.log10(PIXEL_MAX / np.sqrt(mse))


# 每个分块的大小为(8,8)
block_shape = (8, 8)

alpha = 2

# 以灰度图模式读取水印图像
wm = cv2.imread('bupt10.bmp', cv2.IMREAD_GRAYSCALE)

# 从二维矩阵转为一维并二值化
wm = wm.flatten() > 128

# 读取载体图像
img = cv2.imread('lena.bmp')
height, width = img.shape[:2]

# 检查水印数据是否大于载体图像的最大可容纳量
assert (height * width) // (block_shape[0] * block_shape[1]) >= len(wm)

# 取图像的蓝色分量来隐藏
img_b, img_g, img_r = cv2.split(img)

# 对图像进行分块
img_b_blocks = img_to_blocks(img_b, block_shape)

# 信息嵌入
img_b_blocks_embedded = img_b_blocks.copy()
for i in range(len(wm)):
    block = img_b_blocks[i]

    # 对图像分块进行DCT变换
    block_dct = cv2.dct(block)
    block_dct_embedded = block_dct.copy()

    # 选择(4,1)和(3,2)这一对系数
    if wm[i] == 0 and block_dct_embedded[4][1] <= block_dct_embedded[3][2]:
        block_dct_embedded[4][1], block_dct_embedded[3][2] = block_dct_embedded[3][2], block_dct_embedded[4][1]
        # 将原本小的系数调整更小，使得系数差别变大
        block_dct_embedded[3][2] -= alpha

    elif wm[i] == 1 and block_dct_embedded[4][1] >= block_dct_embedded[3][2]:
        block_dct_embedded[4][1], block_dct_embedded[3][2] = block_dct_embedded[3][2], block_dct_embedded[4][1]
        # 将原本小的系数调整更小，使得系数差别变大
        block_dct_embedded[4][1] -= alpha

    # 对图像分块进行DCT逆变换
    block_embedded = cv2.idct(block_dct_embedded)
    img_b_blocks_embedded[i] = block_embedded

# 将分块合并为完整的图像蓝色分量
img_b_embedded = blocks_to_img(img_b_blocks_embedded, img.shape[:2])

# 与图像的绿色、红色分量合并
img_embedded = cv2.merge([img_b_embedded.astype(np.uint8), img_g, img_r])

# 计算峰值信噪比
psnr = PSNR(img, img_embedded)
print(f'PSNR = {round(psnr,3)} dB')

# 保存嵌入水印的载体图像
cv2.imwrite('lenastegoB.bmp', img_embedded)

# 展示原图像和嵌入图像
plt.figure(figsize=(10, 6))

plt.subplot(121)
plt.imshow(cv2.cvtColor(img, cv2.COLOR_BGR2RGB))
plt.title('Original Image 2019212473')
plt.xticks([]), plt.yticks([])

plt.subplot(122)
plt.imshow(cv2.cvtColor(img_embedded, cv2.COLOR_BGR2RGB))
plt.title('Embedded Image 2019212473')
plt.xticks([]), plt.yticks([])

plt.show()

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95

实验结果：
（2）计算嵌入秘密信息后lenastegoB.bmp的峰值信噪比。输出峰值信噪比，要求峰值信噪比精确到小数点后三位。

代码是作业（1）中代码。

实验结果：

（3）在没有任何攻击的情况下从lenastegoB.bmp中图像提取水印二值图像bupt1.bmp，计算bupt10.bmp和bupt1.bmp的归一化函数NC，在没有任何攻击的情况下，提取出来的水印信息和原始的水印信息完全相同，归一化函数的值为1。归一化的计算方法参考Wav音频的LSB信息隐藏。

实验步骤：
（1）首先确定分块大小，设置水印图像的宽高，读取携密图像
（2）取图像的一层来提取，对图像进行分块，并对蓝色分量进行信息提取
（3）信息提取需要进行DCT变换，将其从一维转为二维矩阵，以灰度图模式读取水印图像
（4）计算NC值，保存提取的水印图像，并展示嵌入图像、水印原图像和提取出的水印图像

代码：

# 导入图像处理库
import cv2
# 导入数学计算库
import numpy as np
# 导入绘图库
import matplotlib.pyplot as plt
# 从block.py导入图像分块的函数
from block import img_to_blocks


# 计算NC值的函数
def NC(template, img):
    template = template.astype(np.uint8)
    img = img.astype(np.uint8)
    return cv2.matchTemplate(img, template, cv2.TM_CCORR_NORMED)[0][0]


# 每个分块的大小为(8,8)
block_shape = (8, 8)

# 设置水印图像的宽高
wm_height = 10
wm_width = 10

# 读取携密图像
img = cv2.imread('lenastegoB.bmp')
height, width = img.shape[:2]

# 取图像的一层来提取
img_b, img_g, img_r = cv2.split(img)

# 对图像进行分块
img_b_blocks = img_to_blocks(img_b, block_shape)

# 信息提取
wm = np.zeros(wm_height * wm_width, dtype=np.uint8)
for i in range(len(wm)):
    block = img_b_blocks[i]

    # 对图像分块进行DCT变换
    dct_block = cv2.dct(block)

    if dct_block[4][1] <= dct_block[3][2]:
        wm[i] = 255
    else:
        wm[i] = 0

# 从一维转为二维矩阵
wm = np.reshape(wm, (wm_height, wm_width))

# 以灰度图模式读取水印图像
wm_original = cv2.imread('bupt10.bmp', cv2.IMREAD_GRAYSCALE)

# 计算NC值
nc = NC(wm_original, wm)
print(f'NC = {nc * 100} %')

# 保存提取出的水印图像
cv2.imwrite('bupt1.bmp', wm)

# 展示嵌入图像、水印原图像和提取出的水印图像
plt.figure(figsize=(15, 6))

plt.subplot(131)
plt.imshow(cv2.cvtColor(img, cv2.COLOR_BGR2RGB))
plt.title('Embedded Image 2019212473')
plt.xticks([]), plt.yticks([])

plt.subplot(132)
plt.imshow(wm_original, 'gray')
plt.title('Original Watermark 2019212473')
plt.xticks([]), plt.yticks([])

plt.subplot(133)
plt.imshow(wm, 'gray')
plt.title('Extracted Watermark 2019212473')
plt.xticks([]), plt.yticks([])

plt.show()

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80

实验结果：


（4）在lenastegoB.bmp图像中添加白噪声，将图像另存为lenastegoBl.bmp。从
添加白噪声后的携密图像lenastegoB1.bmp中提取水印二值图像bupt2.bmp，计算bupt10.bmp和bup2.bmp的归一化函数NC，因为携密图像遭受了高斯噪声攻击，提取出来的水印信息和原始的水印信息完全相同，归一化函数的值肯定不为1。计算这个NC值并输出。要求NC的值精确到小数点后三位。

实验步骤：
（1）首先确定分块大小，设置水印图像的宽高，读取携密图像
（2）对携密图像进行高斯噪声攻击，标准差设置为0.002
（3）取图像的一层来提取图像，进行信息提取，需要对图像进行DCT变换
（4）将其从一维变为二维矩阵，以灰度图的模式读取水印图像，并计算NC值
（5）保存提出的水印图像，展示嵌入图像、加入高斯噪声的嵌入图像、水印原图像和提取出的水印图像

代码：

# 导入图像处理库
import cv2
# 导入数学计算库
import numpy as np
# 导入绘图库
import matplotlib.pyplot as plt
# 从block.py导入图像分块的函数
from block import img_to_blocks


# 计算NC值的函数
def NC(template, img):
    template = template.astype(np.uint8)
    img = img.astype(np.uint8)
    return cv2.matchTemplate(img, template, cv2.TM_CCORR_NORMED)[0][0]


# 高斯噪声攻击的函数
def gaussian_attack(img, mean, sigma):
    img = img.astype(np.float32) / 255
    noise = np.random.normal(mean, sigma, img.shape)
    img_gaussian = img + noise
    img_gaussian = np.clip(img_gaussian, 0, 1)
    img_gaussian = np.uint8(img_gaussian * 255)
    return img_gaussian


# 每个分块的大小为(8,8)
block_shape = (8, 8)

# 设置水印图像的宽高
wm_height = 10
wm_width = 10

# 读取携密图像
img = cv2.imread('lenastegoB.bmp')
height, width = img.shape[:2]

# 对携密图像进行高斯噪声攻击，标准差设置为0.002
img_gaussian = gaussian_attack(img, 0, 0.002)

# 取图像的一层来提取
img_b, img_g, img_r = cv2.split(img_gaussian)

# 对图像进行分块
img_b_blocks = img_to_blocks(img_b, block_shape)

# 信息提取
wm = np.zeros(wm_height * wm_width, dtype=np.uint8)
for i in range(len(wm)):
    block = img_b_blocks[i]

    # 对图像分块进行DCT变换
    dct_block = cv2.dct(block)

    if dct_block[4][1] <= dct_block[3][2]:
        wm[i] = 255
    else:
        wm[i] = 0

# 从一维转为二维矩阵
wm = np.reshape(wm, (wm_height, wm_width))

# 以灰度图模式读取水印图像
wm_original = cv2.imread('bupt10.bmp', cv2.IMREAD_GRAYSCALE)

# 计算NC值
nc = np.around(NC(wm_original, wm) * 100, 3)
print(f'NC = {nc} %')

# 保存提取出的水印图像
cv2.imwrite('bupt2.bmp', wm)

# 展示嵌入图像、加入高斯噪声的嵌入图像、水印原图像和提取出的水印图像
plt.figure(figsize=(17, 6))

plt.subplot(141)
plt.imshow(cv2.cvtColor(img, cv2.COLOR_BGR2RGB))
plt.title('Embedded Image')
plt.xticks([]), plt.yticks([])

plt.subplot(142)
plt.imshow(cv2.cvtColor(img_gaussian, cv2.COLOR_BGR2RGB))
plt.title('Gaussian Attacked Embedded Image')
plt.xticks([]), plt.yticks([])

plt.subplot(143)
plt.imshow(wm_original, 'gray')
plt.title('Original Watermark 2019212473')
plt.xticks([]), plt.yticks([])

plt.subplot(144)
plt.imshow(wm, 'gray')
plt.title('Extracted Watermark 2019212473')
plt.xticks([]), plt.yticks([])

plt.show()

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98

实验结果：

五、 实验心得体会
本次实验在进行第一个实验的时候，遇到了block包的报错，才发现未引入block.py代码。
此外，本实验在实验2和实验4尝试了两种不同的方法来保留小数点后三位。实验2中，使用了round()内置函数。round()函数中的第一个变量为需要保留指定位数的小数，第二个变量为指定要保留的位数(该参数默认为0，即不保留小数，直接近似为整数)。实验4中，使用了np.around()函数。np.around()函数对输入浮点数执行5舍6入，5做特殊处理（小数点最后一位为5的舍入为与其值最接近的偶数值）。