信号处理基于小波变换的音频水印嵌入提取matlab源码

Posted 2021-05-18 Matlab走起

大多数自然信号（包括声音和图像）的能量都集中在余弦变换后的低频部分。由于人眼对于细节信息不是很敏感，因此信息含量更少的高频部分可以直接去掉，从而在后续的压缩操作中获得较高的压缩比。

• 1. 一维DCT变换
• 2. 二维DCT变换
• 3. DCT反变换
• 4. 怎样提高二维DCT变换的速度

---

1. 一维DCT变换

一维DCT变换共有8种，最实用的一种公式如下：

其中N是一维数据的元素总数，c(u)系数使得DCT变换矩阵成为正交矩阵，正交特性在二维DCT变换中更能体现其优势。一维DCT变换的复杂度是O(n^2)。

2. 二维DCT变换

二维DCT变换公式如下：

我们将公式变换一下：

可以发现，二维DCT变换其实是在一维DCT变换的基础上，再做一次一维DCT变换。二维DCT也可以写成矩阵相乘的形式：

矩阵A和f必须为长宽均为N的方阵，f不是方阵的情况可以补零再对其做DCT变换。生成矩阵A的matlab代码如下：

二维DCT变换的复杂度达到O(n^4)，所以进行DCT变换的矩阵不宜过大。在实际处理图片的过程中，需要先把矩阵分块，一般分为8*8或16*16大小，这样DCT变换不至于耗费过多的时间。

3. DCT反变换

DCT变换是无损并可逆的，反变换的公式如下：

同样，DCT反变换也可以写成矩阵相乘的形式：

DCT反变换的时间复杂度与DCT变换相同。

3 算 法

本算法采用DWT，包含水印嵌入、水印检测以及水印攻击三个主要部分。水印的工作原理如图2所示。水印的检测需要原始音频数据。



3.1 水印嵌入算法

（1）将待嵌入水印图像置乱。本算法简单地采用伪随机数算法消除数据的相关性。（2）将原始音频数据进行多尺度一维分解，并分别提取低频系数和三层高频系数。为获取较好的稳健性，本算法将水印数据嵌入到音频数据的第三层高频分量上。（3）根据公式Vw（i）=V（i）+（α+e）×W（i）嵌入水印数据。其中V（i）为音频数据位，W（i）为水印数据位，Vw（i）为嵌入水印后的音频数据位，α为水印嵌入强度，e值作为修正，取值10-20。通过试验，发现α值取为0.004时嵌入水印效果较为理想。（4）将嵌入水印数据后的音频进行IDWT变换，即得到包含水印的音频数据。

3.2 水印检测算法

（1）将含有水印的音频数据进行多尺度一维分解，并提取其三层高频分量系数。

（2）检测算法为嵌入算法的逆过程，需要原始音频数据参与检测，表达为：

W（i）=（Vw（i）-V（i））/（α+e），其中α以及e的值同嵌入算法中确定的值一致。

（3）步骤（2）得到的W（i）即为提取的一维水印信息序列，将其进行升维处理，可以得到二维图像形式。此结果便是检测输出的水印。

嵌入水印的原图及一次嵌入后提取的水印分别如图3和图4所示。



4 部分试验

为测试本水印系统的性能，对加水印的音频数据进行各类攻击，这里给出部分实验结果。

定义：



Nc作为衡量所提取水印图像与原始水印图像的相似程度，从未被攻击的含水印音频中直接提取的水印与原水印图像相似度高达0.9998。

（1）二选一迫选实验。对事先不知道精确原始音频信号的测试者分别播放原始音频和嵌入水印后的音频，要求测试者指认原始音频。根据L.BONey等[5]的结论，如果二类音频被指认为原始音频的比例大致相当，则可认为水印嵌入后没有引起人耳感知上的显着差别。试验中随机选取同实验室学生8人，通过分别在不同wav文件中嵌入水印并随机询问的方法，结果约有53.4%的接受询问者认为原始音频音质更好。说明通过本系统嵌入的水印没有引起原始音频音质上的显着改变。

（2）将音频截掉全部数据的n/10（n=1，2，3……），原始音频数据位长度稍大于40 000，从第20 000位开始剪切。

根据图5、图6、图7可知将音频剪切掉约三分之一内容后，仍可提取出较为明显的水印图案。如果剪切部分再多一些，则无法满意地检测出水印。但由于三分之一的剪切率将同时导致载体音频数据的大量丢失，故这个结果是可以接受的。



（3）MP3压缩。目前对音频信号进行MP3压缩编码是较为常用的一种音频处理技术，其目标为在不影响原始音频信号品质的前提下尽可能地减少音频数据量。不同的比特率对应了不同的MP3压缩比。本试验对上面含有水印的一段音频先进行码率为96Kbps的压缩（压缩比为7.4：1），然后进行相映解码处理，检测得到的水印图像如图8所示。



5 结 论

近年来音频数字水印领域尤其是变换域音频水印嵌入与检测方面的研究工作发展迅速，而离散小波分析（DWT）是近年来整个数字水印系统研究的热点之一。算法经实验表明具有很好的隐蔽性，对原始音频的质量几乎没有削弱，具备一定的抵抗剪切攻击及其他攻击的能力。为进一步提高算法的稳健性，应该进一步考虑如何利用更多HAS特性以及水印嵌入的位置和强度。顾及算法实用性，应当考虑增加嵌入水印的容量问题。这些都是需进一步改善的方向。

clear all;
clc;
key=35;
%Arnold置换次数，作为密钥
Orignalmark=double(imread('suda64.bmp'));  %读入64*64的水印图片
[wrow,wcol]=size(Orignalmark);  
if wrow~=wcol 
    error('wrow~=wcol error');
end
%--- 测试密钥key是否超出范围---------
n=check_arnold(wrow);
if (key+1)>n
    error('arnold key error');
end
Arnoldw=arnold(Orignalmark,wrow,key); %对水印图像进行Arnold转化
[X,fs,bits]=wavread('laile.wav'); %读入音频文件
%X=imread('lena.bmp');
figure;    
subplot(2,1,1); 
%imshow(X),title('原始图像');
%X=double(X);
plot(X);     %显示音频文件波形
title('原始音频信号'); 
%sound(X,fs,bits);
%pause;
%水印嵌入--------------------------------------------------
[c,l]=wavedec(X,2,'db4'); %用db4小波对读入的声音文件进行2级小波分解 
ca2=appcoef(c,l,'db4',2); %提取2级小波分解的低频系数和高频系数 
cd2=detcoef(c,l,2); 
cd1=detcoef(c,l,1); 
lca=length(ca2);  %低频长度
blocksize=fix(lca/(wrow*wcol)); %每块的大小
water_vector=reshape(Arnoldw,1,wrow*wcol);  %将置乱后的水印转化为一维的
wlength=wrow*wcol;  %水印的长度
a=0.25;  %量化步长
j=1;
for i=1:wlength
    Block=ca2(j:j+blocksize-1);
    [U,S,V]=svd(double(Block)); 
    cc=floor(S(1,1)/a);  
    if(Arnoldw(i)==1)           %嵌入奇数倍
        if(mod(cc,2)==0)
            cc=cc+1;
        end
        S(1,1)=a*cc;
    end
    if(Arnoldw(i)==0)            %嵌入偶数倍
        if(mod(cc,2)==1)  
            cc=cc+1;
        end
        S(1,1)=a*cc;
    end
    Blockw=U*S*V';          %SVD 逆变换还原 
    ca2(j:j+blocksize-1)=Blockw; 
    j=j+blocksize;
end
c1=[ca2',cd2',cd1']';
 
 
 

完整代码添加QQ1575304183