兩句閑話
老師在課上講了許多圖片隱寫和隱寫分析的方法，在這里我整合一下，並對部分進行代碼實現。
LSB隱寫
LSB隱寫是最基礎、最簡單的隱寫方法，具有容量大、嵌入速度快、對載體圖像質量影響小的特點。
LSB的大意就是最低比特位隱寫。我們將深度為8的BMP圖像，分為8個二值平面（位平面），我們將待嵌入的信息（info）直接寫到最低的位平面上。換句話說，如果秘密信息與最低比特位相同，則不改動；如果秘密信息與最低比特位不同，則使用秘密信息值代替最低比特位。
具體實現如下
from PIL import Image
import math
class LSB:
def __init__(self):
self.im=None
def load_bmp(self,bmp_file):
self.im=Image.open(bmp_file)
self.w,self.h=self.im.size
self.available_info_len=self.w*self.h # 不是絕對可靠的
print ("Load>> 可嵌入",self.available_info_len,"bits的信息")
def write(self,info):
"""先嵌入信息的長度，然后嵌入信息"""
info=self._set_info_len(info)
info_len=len(info)
info_index=0
im_index=0
while True:
if info_index>=info_len:
break
data=info[info_index]
x,y=self._get_xy(im_index)
self._write(x,y,data)
info_index+=1
im_index+=1
def save(self,filename):
self.im.save(filename)
def read(self):
"""先讀出信息的長度，然后讀出信息"""
_len,im_index=self._get_info_len()
info=[]
for i in range(im_index,im_index+_len):
x,y=self._get_xy(i)
data=self._read(x,y)
info.append(data)
return info
#===============================================================#
def _get_xy(self,l):
return l%self.w,int(l/self.w)
def _set_info_len(self,info):
l=int(math.log(self.available_info_len,2))+1
info_len=[0]*l
_len=len(info)
info_len[-len(bin(_len))+2:]=[int(i) for i in bin(_len)[2:]]
return info_len+info
def _get_info_len(self):
l=int(math.log(self.w*self.h,2))+1
len_list=[]
for i in range(l):
x,y=self._get_xy(i)
_d=self._read(x,y)
len_list.append(str(_d))
_len=''.join(len_list)
_len=int(_len,2)
return _len,l
def _write(self,x,y,data):
origin=self.im.getpixel((x,y))
lower_bit=origin%2
if lower_bit==data:
pass
elif (lower_bit,data) == (0,1):
self.im.putpixel((x,y),origin+1)
elif (lower_bit,data) == (1,0):
self.im.putpixel((x,y),origin-1)
def _read(self,x,y):
data=self.im.getpixel((x,y))
return data%2
if __name__=="__main__":
lsb=LSB()
# 寫
lsb.load_bmp('test.bmp')
info1=[0,1,0,1,1,0,1,0]
lsb.write(info1)
lsb.save('lsb.bmp')
# 讀
lsb.load_bmp('lsb.bmp')
info2=lsb.read()
print (info2)
在這里，我們定義幾個指標來評價隱寫算法。
假設，某灰度圖像的大小為$M\times N$，深度為8。
嵌入容量 $M\times N$bit。
嵌入率 $\frac{嵌入容量}{總容量}$，LSB的嵌入率為$\frac{1}{8}=12.5\%$。
MSE mean square error，平均方根誤差，$MSE=\frac{\sum^{M}_{m=1}\sum^{N}_{n-1}d(m,n)^2}{M \times N}$，這里的$d(m,n)$指的是，原圖像和修改后的圖像在$(m,n)$位置上的像素點之差。
PSNR peak signal-to-noise ratio，峰值信噪比，$PSNR=-10log\{\frac{MSE}{255^2MN}\}$。
更多的評價方法還有，VQM(vedio quality measurement), SSIM(structural similarity index)等。
面向JPEG的圖像隱寫（1）：Jsteg隱寫
關於JPEG格式，可以看看這篇博客 JPEG圖像壓縮算法流程詳解。JPEG壓縮中，最主要的就是DCT變換。
Jsteg隱寫是將秘密信息嵌入在量化后的DCT系數的LSB上，但原始值為-1,0，+1的DCT系數除外。此外，由於量化后的DCT系數中有負數，編程的時候需要格外注意以下。
具體實現如下
import math
class Jsteg:
def __init__(self):
self.sequence_after_dct=None
def set_sequence_after_dct(self,sequence_after_dct):
self.sequence_after_dct=sequence_after_dct
self.available_info_len=len([i for i in self.sequence_after_dct if i not in (-1,1,0)]) # 不是絕對可靠的
print ("Load>> 可嵌入",self.available_info_len,'bits')
def get_sequence_after_dct(self):
return self.sequence_after_dct
def write(self,info):
"""先嵌入信息的長度，然后嵌入信息"""
info=self._set_info_len(info)
info_len=len(info)
info_index=0
im_index=0
while True:
if info_index>=info_len:
break
data=info[info_index]
if self._write(im_index,data):
info_index+=1
im_index+=1
def read(self):
"""先讀出信息的長度，然后讀出信息"""
_len,sequence_index=self._get_info_len()
info=[]
info_index=0
while True:
if info_index>=_len:
break
data=self._read(sequence_index)
if data!=None:
info.append(data)
info_index+=1
sequence_index+=1
return info
#===============================================================#
def _set_info_len(self,info):
l=int(math.log(self.available_info_len,2))+1
info_len=[0]*l
_len=len(info)
info_len[-len(bin(_len))+2:]=[int(i) for i in bin(_len)[2:]]
return info_len+info
def _get_info_len(self):
l=int(math.log(self.available_info_len,2))+1
len_list=[]
_l_index=0
_seq_index=0
while True:
if _l_index>=l:
break
_d=self._read(_seq_index)
if _d!=None:
len_list.append(str(_d))
_l_index+=1
_seq_index+=1
_len=''.join(len_list)
_len=int(_len,2)
return _len,_seq_index
def _write(self,index,data):
origin=self.sequence_after_dct[index]
if origin in (-1,1,0):
return False
lower_bit=origin%2
if lower_bit==data:
pass
elif origin>0:
if (lower_bit,data) == (0,1):
self.sequence_after_dct[index]=origin+1
elif (lower_bit,data) == (1,0):
self.sequence_after_dct[index]=origin-1
elif origin<0:
if (lower_bit,data) == (0,1):
self.sequence_after_dct[index]=origin-1
elif (lower_bit,data) == (1,0):
self.sequence_after_dct[index]=origin+1
return True
def _read(self,index):
if self.sequence_after_dct[index] not in (-1,1,0):
return self.sequence_after_dct[index]%2
else:
return None
if __name__=="__main__":
jsteg=Jsteg()
# 寫
sequence_after_dct=[-1,0,1]*100+[i for i in range(-7,500)]
jsteg.set_sequence_after_dct(sequence_after_dct)
info1=[0,1,0,1,1,0,1,0]
jsteg.write(info1)
sequence_after_dct2=jsteg.get_sequence_after_dct()
# 讀
jsteg.set_sequence_after_dct(sequence_after_dct2)
info2=jsteg.read()
print (info2)
在上面，我們實現了對量化后的DCT系數的隱寫。至於如何得到DCT系數，可以使用opencv中的函數，如下
import cv2
import numpy as np
def dct(m):
m = np.float32(m)/255.0
return cv2.dct(m)*255
面向JPEG的圖像隱寫（2）：F3隱寫
在Jsetg隱寫方法中，原始值為-1，0，+1的DCT系數，不負載秘密信息，但是量化后的DCT系數中卻有大量的-1，0，+1（以0居多），這說明Jsetg的嵌入率會很小。為了改善這一狀況，人們提出了F3隱寫。
F3則對原始值為+1和-1的DCT系數，進行了利用。F3隱寫的規則如下
- (1) 每個非0的DCT數據用於隱藏1比特秘密信息，為0的DCT系數不負載秘密信息。
- (2) 如果秘密信息與DCT的LSB相同，便不作改動；如果不同，將DCT系數的絕對值減小1，符號不變。
- (3) 當原始值為+1或-1且預嵌入秘密信息為0時，將這個位置歸0並視為無效，在下一個DCT系數上重新嵌入。
我們可以看出來，F3對Jsteg的改動並不大。因此，在代碼實現上，我們可以復用Jsteg的代碼，具體如下
from jsteg import Jsteg
import math
class F3(Jsteg):
def __init__(self):
Jsteg.__init__(self)
def set_sequence_after_dct(self,sequence_after_dct):
self.sequence_after_dct=sequence_after_dct
sum_len=len(self.sequence_after_dct)
zero_len=len([i for i in self.sequence_after_dct if i==0])
one_len=len([i for i in self.sequence_after_dct if i in (-1,1)])
self.available_info_len=sum_len-zero_len-one_len # 不是特別可靠
print ("Load>> 大約可嵌入",sum_len-zero_len-int(one_len/2),'bits')
print ("Load>> 最少可嵌入",self.available_info_len,'bits\n')
def _write(self,index,data):
origin=self.sequence_after_dct[index]
if origin == 0:
return False
elif origin in (-1,1) and data==0:
self.sequence_after_dct[index]=0
return False
lower_bit=origin%2
if lower_bit==data:
pass
elif origin>0:
self.sequence_after_dct[index]=origin-1
elif origin<0:
self.sequence_after_dct[index]=origin+1
return True
def _read(self,index):
if self.sequence_after_dct[index] != 0:
return self.sequence_after_dct[index]%2
else:
return None
if __name__=="__main__":
f3=F3()
# 寫
sequence_after_dct=[-1,0,1]*100+[i for i in range(-7,500)]
f3.set_sequence_after_dct(sequence_after_dct)
info1=[0,1,0,1,1,0,1,0]
f3.write(info1)
sequence_after_dct2=f3.get_sequence_after_dct()
# 讀
f3.set_sequence_after_dct(sequence_after_dct2)
info2=f3.read()
print (info2)
可嵌入容量的計算
這里需要說的是，由於F3隱寫特殊的規則，我們無法精確得到可嵌入的信息的容量，我們只能得到最小值，即原始值為非0，-1，+1的像素點的數量。但是，我們可以得到一個數學期望。但是這個期望等於多少呢？我們來算一下。
為了嚴謹性，我們先列出幾條假設：
- (1) 待嵌入信息為01串。在此01串中，0和1隨機均勻分布，且0和1出現的概率分別為50%。
- (2) 假設系數表中不同系數的出現是隨機的，我們忽略它們出現的次序，如非0、-1、+1的出現總是相鄰的。
此外，我們設量化后的DCT系數表中，0的概率為$p_0$，-1和1的概率為$p_1$，其他數字出現的概率為$p_2$；DCT系數表的長度為n；待嵌入的01串的長度為m。
則假設我們要嵌入0，我們需要DCT系數的個數為$\frac{1}{p_2}$；假設我們要嵌入1，我們需要DCT系數的個數是$\frac{1}{p_1+p_2}$。由於01出現的概率分別為50%，因此我們嵌入一位所需要的DCT系數的個數為$\frac{\frac{1}{p_2}+\frac{1}{p_1+p_2}}{2}=\frac{p_1+2p_2}{2p_2(p_1+p_2)}$，因此m與n的關系為$m=\frac{n}{\frac{p_1+2p_2}{2p_2(p_1+p_2)}}=\frac{2np_2(p_1+p_2)}{p_1+2p_2}$。
在實際實施的，我們統計0的數量$n_0$，-1和1的數量$n_1$，其他數字出現的概率為$n_2$，於是$m=\frac{2n_2(n_1+n_2)}{n_1+2n_2}$。
面向JPEG的圖像隱寫（3）：F5隱寫
調色板隱寫(EZStego隱寫)
首先，介紹調色板圖像。
調色板圖像 調色板圖像是互聯網上常見的一種圖像格式，其中含有一個不超過256種顏色的調色板，並定義了每種顏色對應的R,G,B各顏色分量值，圖像內容中的每個像素是不超過8比特信息的一個索引值，其指向的調色板中的對應顏色即該像素中的真實顏色。常見的調色板圖像格式是GIF,PNG。
EZStego隱寫
- (1) 將調色板的顏色亮度依次排序，其中顏色的亮度由不同的顏色分量線性疊加而成，其表達式為$Y=0.299R+0.587G+0.114B$。
- (2) 為每個顏色分配一個亮度序號。
- (3) 將調色板圖像像素內容使用LSB隱寫代替，並將圖像像素索引值改為新的亮度序號所對用的索引值。
- (4) 用奇數序號表示嵌入秘密比特1，用偶數序號表示嵌入秘密比特0。
python代碼
from PIL import Image
import numpy as np
import math
"""
假設調色板索引為
0 1 2 3 4 5 6 7
假設亮度序號(Y_index)為
index: 0 1 2 3 4 5 6 7
2 5 4 1 7 3 6 0
則
# Y_index_inverse
index: 0 1 2 3 4 5 6 7
7 3 0 5 2 1 6 4
# 例子
載體 [3 0 6 4] ; 待嵌入信息 0110
嵌入：[3 0 6 4]=>[5 7 6 2]=>(by 0110) [4 7 7 2]=>[7 0 0 4]
結果：[3 0 6 4]=>(by 0110)=>[7 0 0 4]
提取：[7 0 0 4]=>[4 7 7 2]=>0110
"""
class GIF_Steg:
def __init__(self):
self.im=None
def load_gif(self,gif_file):
self.im=Image.open(gif_file)
self._load_palette()
self._sort_palette()
self._load_palette_data()
self.available_info_len=len(self.palette_data)
def write(self,info):
info=self._set_info_len(info)
self.palette_data=self._write(self.palette_data,info)
def read(self):
_len,im_index=self._get_info_len()
info=self._read(self.palette_data[im_index:im_index+_len])
return info
def save(self,filename):
self.im.save(filename)
#==========================================#
def _load_palette(self):
self.palette=[]
palette=self.im.palette.palette
for i in range(int(len(palette)/3)):
self.palette.append((palette[3*i],palette[3*i+1],palette[3*i+2]))
def _sort_palette(self):
f=lambda t:0.299*t[0]+0.587*t[1]+0.114*t[2]
Y=[f(t) for t in self.palette]
self.Y_index=np.argsort(Y)
self.Y_index_inverse=[0]*256
for i in range(len(self.Y_index)):
self.Y_index_inverse[self.Y_index[i]]=i
def _load_palette_data(self):
self.palette_data=self.im.getpalette()
def _set_info_len(self,info):
l=int(math.log(self.available_info_len,2))+1
info_len=[0]*l
_len=len(info)
info_len[-len(bin(_len))+2:]=[int(i) for i in bin(_len)[2:]]
return info_len+info
def _get_info_len(self):
l=int(math.log(self.available_info_len,2))+1
len_list=[]
for i in range(l):
_d=self._get_lsb(self.palette_data[i])
len_list.append(str(_d))
_len=''.join(len_list)
_len=int(_len,2)
return _len,l
def _write(self,palette_data,info):
for i in range(len(info)):
Y_index=self.Y_index_inverse[palette_data[i]]
lower_bit=Y_index%2
if lower_bit==info[i]:
pass
elif (lower_bit,info[i])==(0,1):
palette_data[i]=self.Y_index[Y_index+1]
elif (lower_bit,info[i])==(1,0):
palette_data[i]=self.Y_index[Y_index-1]
return palette_data
def _read(self,palette_data):
info=[]
for i in range(len(palette_data)):
info.append(self._get_lsb(palette_data[i]))
return info
def _get_lsb(self,_palette_data):
return self.Y_index_inverse[_palette_data]%2
if __name__=="__main__":
gs=GIF_Steg()
gs.load_gif('4.1.05.gif')
gs.write([0,1,1,0,0,1,0,0,0,0,0])
print (gs.read())
BPCS隱寫
PVD隱寫
卡方分析
RS分析
RS隱寫分析的原理，在百度文庫的這個ppt上說的比較清楚。
首先介紹像素翻轉$F_1,F_0,F_{-1}$。$F_1$是像素值$2n$與$2n+1$之間的變換，$F_{-1}$是像素值$2n-1$與$2n$之間的變換，$F_0$則是像素值不發生改變。即
設一掩碼算子$m=(m_1,m_2,\cdots,m_n),(m_i\in {0,1})$。現在定義$F_m$與$F_{-m}$。對於長度為n的像素值的序列$G$，$F_m(G)=(F_{m_1}(G[1]),\cdots,F_{m_i}(G[i]),\cdots,F_{m_n}(G[n]))$。相應地，$F_{-m}(G)=(F_{-m_1}(G[1]),\cdots,F_{-m_i}(G[i]),\cdots,F_{-m_n}(G[n]))$。
現在，我們定義像素相關性，設$G$長度為n的像素值的序列，$G=(x_1,x_2,\cdots,x_n)$。則$序列G$像素相關性$f(G)=\sum_{i=1}^{n-1}{|x_{i+1}-x_i|}$。
大量實驗表明，當一個像素值序列經歷$F_m$或$F_{-m}$之后，像素相關性的變化會隨着圖片中嵌入秘密信息的數量會呈現出一些規律。
我們將圖片分塊，每一塊通過Z字形掃描變成一段序列，這樣我們就得到了多個像素點序列。對所有序列使用非負翻轉$F_m$和$F_{-m}$翻轉，像素相關性增加或減少的比例，我們分別設為$R_m,S_m,R_{-m},S_{-m}$。即
- $R_m$ 為$F_m$作用下像素相關性增加占所有像素組的比例
- $R_{-m}$ 為$F_{-m}$作用下像素相關性增加占所有像素組的比例
- $S_m$ 為$F_{m}$作用下像素相關性減少占所有像素組的比例
- $S_{-m}$ 為$F_{-m}$作用下像素相關性減少占所有像素組的比例
假設一圖像嵌入了秘密信息，嵌入率為$\alpha$，即原圖中比例為$\frac{\alpha}{2}$的像素值發生了改變，那么$\alpha$與$R_m,R_{-m},S_m,S_{-m}$的關系如下圖(大量實驗的結果)
下面的代碼能夠得到上圖。至於如何根據一張隱寫的圖片，得到圖片中是否經過隱寫，以及得到嵌入率，這就是另一個問題。這個問題比較復雜，上面給的ppt中，有介紹。
import sys
import math
import numpy as np
from PIL import Image
import random
def get_index_matrix(n):
"""
得到n階zigzag掃描矩陣
"""
I=np.array(range(n))
J=I.reshape(-1,n).T
M=((I+J)*(I+J+2)+(I-J)*(-1)**(I+J))/2
one_tril=np.triu(np.ones((n,n)))[:,::-1]
M=M*one_tril
M=M+(n**2-1-M)[::-1,::-1]*(1-one_tril)
return M.astype(int)
def get_mask(n):
"""
得到掩碼m
"""
return np.random.randint(low=0,high=2,size=n)
class RS:
def __init__(self):
self._region_length=8
self.set_parameter()
def load_bmp(self,bmp_file):
"""
加載bmp文件
"""
self.im=Image.open(bmp_file)
self.w,self.h=self.im.size
print (">> 加載圖片，圖片尺寸：",self.w,"x",self.h)
def set_parameter(self,_region_length=8):
self._region_length=_region_length
self._zigzag_index_matrix=get_index_matrix(_region_length)
self._m=get_mask(_region_length**2)
def analyse(self):
_rs1=[0,0,0,0] # [Rm,Sm,R-m,S-m]
_rs2=[0,0,0,0] # [Rm,Sm,R-m,S-m]
self._RS_build(_rs1,_rs2)
_sum=math.ceil(self.w/self._region_length)*math.ceil(self.h/self._region_length)
_rs1=[i/_sum for i in _rs1]
_rs2=[i/_sum for i in _rs2]
res=self._get_insert_rate(_rs1,_rs2)
print (res)
############## unfinished
def get_RS_map(self,n=100):
"""
得到點集 (嵌入率-RS)
"""
res=[]
for i in range(n+1):
_rs=[0,0,0,0]
rate=i/n
self._RS_build_by_rate(_rs,rate)
print (rate,_rs)
res.append((rate,_rs))
return res
######################################################
def _RS_build(self,_rs1,_rs2):
row=math.ceil(self.w/self._region_length)
column=math.ceil(self.h/self._region_length)
for i in range(row):
for j in range(column):
# 從圖像取出一塊區域，進行zigzag掃描
box=np.array([i,j,i+1,j+1])*self._region_length
region=self.im.crop(box)
region=np.array(region)
sequence=self._zigzagScan(region)
# 對RS進行統計
self._rs_build(sequence,_rs1)
# 進行正翻轉，得到修改率為1-a/2的序列，對RS進行統計
sequence=self._Fm(sequence,np.ones(self._region_length**2).astype(int))
self._rs_build(sequence,_rs2)
def _RS_build_by_rate(self,_rs,rate):
"""
根據嵌入率得到RS的值
"""
row=math.ceil(self.w/self._region_length)
column=math.ceil(self.h/self._region_length)
for i in range(row):
for j in range(column):
# 從圖像取出一塊區域，進行zigzag掃描
box=np.array([i,j,i+1,j+1])*self._region_length
region=self.im.crop(box)
region=np.array(region)
sequence=self._zigzagScan(region)
# 以概率rate,嵌入01
sequence=self._random_inject(sequence,rate)
# 對RS進行統計
self._rs_build(sequence,_rs)
def _zigzagScan(self,m):
"""
Z字形掃描
"""
sequence = np.zeros(self._region_length**2,).astype(int)
for i in range(self._region_length):
for j in range(self._region_length):
index = self._zigzag_index_matrix[i][j]
sequence[index] = m[i,j]
return sequence
def _random_inject(self,sequence,rate):
"""
隨機嵌入秘密信息，嵌入率rate
"""
m=np.ceil(np.random.random(self._region_length**2)-rate/2).astype(int)
return self._Fm(sequence,m)
def _rs_build(self,sequence,_rs):
"""
根據sequence修改RS的值
"""
r1=self._get_relativity(sequence)
r2=self._get_relativity(self._Fm(sequence, self._m))
r3=self._get_relativity(self._Fm(sequence,-self._m))
if r1<r2:
_rs[0]+=1
elif r1>r2:
_rs[1]+=1
if r1<r3:
_rs[2]+=1
elif r1>r3:
_rs[3]+=1
def _get_relativity(self,sequence):
"""
得到像素相關性
"""
a=np.abs(np.array(sequence)[1:]-np.array(sequence)[:1])
return np.sum(a)
def _Fm(self,sequence,m):
"""
由m定義的翻轉
"""
# [0,1,-1]
# (x+0)^0-0,(x+0)^1-0,(x-1)^1+1
# ((x+a)^b)-a
a=np.floor(m/2).astype(int)
b=np.abs(m).astype(int)
return ((sequence+a)^b)-a
if __name__=="__main__":
rs=RS()
rs.load_bmp("../_data/misc/5.3.01.tiff")
# rs.analyse()
res=rs.get_RS_map()
import matplotlib.pyplot as plt
for i in range(4):
plt.plot([p[0] for p in res],[p[1][i] for p in res],'ro')
plt.show()
對於上面的代碼，有三點需要說明。
1.如何實現zigzag掃描
在代碼中，我們是通過一個索引矩陣來實現zigzag掃描的，其中，八階的索引矩陣如下
那么，對於其他階數的索引矩陣，怎么得出呢？設M為n階索引矩陣，則有如下的關系
那這個是怎么得到的呢？當然是自己推啦。這里給出推導思路，首先求出$M[0,j]$和$M[i,0]$的關系式，然后利用下面兩個關系式得到$M[i,j]$的表達式
至於有沒有必要，費勁波折得到這個，那我就不知道了。
2.如何快速實現翻轉
使用異或，我們能實現快速的翻轉。$F_0(x)=x\oplus 0,\\F_1(x)=x\oplus 1,\\F_{-1}(x)=((x-1)\oplus 1)+1$。\\為了快速地實現$F_m$，我們定義$F(a,b,x)=((x-a)\oplus b)+a$。於是$F_0(x)=F(0,0,x),F_1(x)=F(0,1,x),F_{-1}(x)=F(1,1,x)$。進一步，對於翻轉$F_i$，令$a=\lfloor i/2\rfloor,b=|i|$,於是$F_i(x)=F(a,b,x)$。
這樣有什么用呢？答案就是大大方便了矩陣運算。不過，應該有比這更快的計算方法。這里不做研究。
3.敏感性分析
m怎么確定？按照n*n分塊，n怎么確定？
在代碼中，我是隨機的生成一個含0和1比例各50%的一個向量。為什么是50%呢？我發現對於有些圖像設為50%會得到比較好的結果(圖像比較合乎規律)，有些圖像比例應該設為90%才會得到比較好的結果。
n的大小又會怎么影響結果？
這里的水，就比較深了。