根据文献中描述的人类视觉系统模型编写数字水印的嵌入程序。可是由于像素临域的计算繁多，加上第三个函数后计算很慢效果还很差，期待高手指点，具体算法在文献中，文献本帖下面的文件下载。
我做的程序有：
%在密钥控制下生成和图象Ａ大小相同的二值｛+1,-1｝二维伪随机图像
function xtheta=pseudorandom(a,k)
[row,col]=size(a);
rand(’state’,k);
temp=randperm(row*col);
b=zeros(row,col);
for i=1:row
　　for j=1:col
　　if temp((i-1)*col+j)>=(row*col/2)
　　b(i,j)=1;
　　else
　　b(i,j)=-1;
　　end
　　end
end

xtheta{1,1}=b(1:round(row/2),round(col/2)+1:col);
xtheta{1,2}=b(round(row/2)+1:row,round(col/2)+1:col);
xtheta{1,3}=b(round(row/2)+1:row,1:round(col/2));
xtheta{1,4}=b(1:round(row/2),1:round(col/2));

%根据人类视觉系统特点，对图像的小波变换系数进行调制
function II0=hvsmodulation(I,I0,I33,xtheta)
[row,col]=size(I0{1,2});
alpha=0.5;

for theta=1:3
　　for i=1:row
　　for j=1:col
　　II0{1,theta}(i,j)=I0{1,theta}(i,j)+alpha*weight(0,theta,i,j,I,I33)*xtheta{1,theta}(i,j);
　　end
　　end
end
%权重函数，
function y=weight(l,theta,i,j,I,I33)

q0=fun1(0,theta)*fun2(0,i,j,I33);%*fun3(0,i,j,I,I33);

y=q0/2;
%权重子函数（１），人眼对高频波段中的噪声不敏感，对对角线方向也是如此
function y=fun1(l,theta)

if theta==2
　　y=sqrt(2);
else
　　y=1;
end

if l==0
　　;
elseif l==1
　　y=y*0.32;
elseif l==2
　　y=y*0.16;
elseif l==3
　　y=y*0.10;
end

%权重子函数（２），人眼对图像中高亮或低亮区域中的噪声较不敏感
function y=fun2(l,i,j,I33)

i=1+floor((i-1)/(2^(3-l)));
j=1+floor((j-1)/(2^(3-l)));
L=I33(i,j);
L=L/256;

if L<0.5
　　L=1-L;
end

y=1+L;
　　
%权重子函数（３），人眼对稠密纹理中的噪声不敏感，但对纹理边缘较敏感
function f=fun3(l,i,j,I,I33)

xv=[0,1];
yv=[0,1];
iv=1+xv+floor((i-1)/(2^(3-l)));
jv=1+yv+floor((j-1)/(2^(3-l)));
if iv(1,2)>16
　　iv(1,:)=16;
end
if jv(1,2)>16
　　jv(1,:)=16;
end
Iv0=I33(jv,iv);
Iv=[Iv0(1,1),Iv0(1,2),Iv0(2,1),Iv0(2,2)];
f3=var(Iv);

f=0;

for k=0:3-l
　　 f1=1/(16^k);
　　 for theta=0:2
　　 for x=0:1
　　 for y=0:1
　　 Ik=I{1,k+l+1};
　　 Iktheta=Ik{1,theta+1};
　　 ii=y+ceil((i)/(2^k));jj=x+ceil((j)/(2^k));
　　 if ii>128/(2^k)|jj>128/(2^k)
　　 f2=0;%Iktheta(128/(2^k),128/(2^k));
　　 else
　　 f2=Iktheta(ii,jj);
　　 end
　　 f2=f2^2;
　　 f=f+f1*f2;
　　 end
　　end
　　end
end
f=(f^0.2)*f3;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
以下是直接在命令区域输入的程序。’lena.gif’可在http://telesun.insa-lyon.fr/~telesun/Img/lena.gif下载。

A=imread(’lena.gif’);X=double(A);
[I03,I02,I00,I01] = dwt2(X,’haar’);
I0={I00,I01,I02,I03};

[I13,I12,I10,I11] = dwt2(I03,’haar’);
I1={I10,I11,I12,I13};

[I23,I22,I20,I21] = dwt2(I13,’haar’);
I2={I20,I21,I22,I23};

[I33,I32,I30,I31] = dwt2(I23,’haar’);
I3={I30,I31,I32,I33};

I={I0,I1,I2,I3};

xtheta=pseudorandom(X,666);
II0=hvsmodulation(I,I0,I33,xtheta);

watermarked_image=idwt2(I03,II0{1,3},II0{1,2},II0{1,1},’haar’);
difference=abs(double(A)-double(watermarked_image));

figure,
subplot(2,2,1);imshow(A,[]);
subplot(2,2,2);imshow(watermarked_image,[]);
subplot(2,2,3);imshow(difference,[]);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
如果只用fun1、fun2，效果还算可以，问题应该出在fun３上。
希望XDJM们耐心看一下，如有指教，不尽感激~~
一个人做，苦于没有交流……

[此贴子已经被作者于2006-5-1 19:23:17编辑过]