1、偏微分方程与图像处理(GAC 的水平集方法)实验二 GAC 的水平集方法一 实验目的采用 GAC 模型的水平集方法检测图像中对象的轮廓,以便有效地进行分割。二 原理分析推广 GAC 模型的水平集方法对应的 PDE 为:(3.31)ugcugkt按照上式,曲线运动将受两种“力”的支配,第一种力来自于曲率几何形变曲率运动() ,不过它的强弱还要受到因子 的影响。k ()gI为图象 I(x,y)的梯度模值,函数 g(r) 是可以是任何具有单调减性的函数。I因为图象梯度模值 在图象的边缘附近有较大值,从而使 g( )取极小的值,故在图II象边缘附近,该作用力将会变的很小,因此有时将边缘函数 称之为边缘
2、停止函数。(g常数 c 的作用是加速曲线向内部收缩。第二种力来自于 g 的梯度 ,它是一种不论当前 C 的局部是在对象内部或外(1,2)部,都能将曲线引向边界的“吸引力” 。从而 总是使曲线向着更接近于边界线的gu方向运动,最终达到贴近对象边界的稳定状态。由于这两种作用使曲线演化可最终达到紧靠轮廓这一稳定状态而不再继续演化。采用单边迎风方案,根据(1.76)式的数值方案实现上式:考虑到 ,0gc可得: (1)()()nnijijijut() () () ()max,0in(1,0max(2,0in(2,0xij xij yij yijDuDuDuDu(2.1)(0)2()21nijxijyij
3、gk其中() 2222a(,)(in,0)(ax,0)(in,0)xij xij yij yijuuuu(2.3)中心差分 (2.2),1,(0)2ijijxijD向前单边差分 (2.3),1,xijijijDu向后单边差分 (2.4),ijijij三 编程过程1 准备工作1)读入图像 I,将其转化为灰度图象,重新调整图象的大小为100,100。2)进行预平滑,即对图象进行高斯卷积滤波。3)计算图像梯度模值 , ,代入函数 计算 ,然后计算 gIrIexp()grk()gI的梯度 。(1,2)g4)选取初始封闭曲线 ,使其从外部全部包住对象,简单起见这里将 选取为以图像中0C0C心为圆心的封闭
4、圆。5)根据 初始化水平集 。00u2 迭代运算1)根据(2.1)式进行迭代运算,由 计算 。()nij(1)niju2)每迭代 5 次,进行一次重新初始化,以免累计误差。具体的方法是根据当前演化得到的检测零水平集则为当前 ,根据当前 重新初始化水平集 。uC3)每 10 次观察一次零水平集,当演化曲线迭代 400 次,则停止迭代。3 参数取值:c=3 4, =0.05 0.1, N=3 5, k=1 2.t四 程序1、主程序:% 用 GAC 水平集演化方法检测图像轮廓close all;clear all;a=imread(3.bmp);% figure;imshow(a);a=rgb2gr
5、ay(a); im=imresize(a,100 100);% figure;imshow(im);imsize=size(im);im_1=double(im);% 对图像进行高斯滤波sigma=2;gauss_filter=fspecial(gaussian,3,sigma); %默认值3*3,SIGMA=0.5b=imfilter(im_1,gauss_filter,conv);% 计算图像梯度Ix,Iy和梯度模值 deltIIx,Iy=gradient(b);deltI=abs(sqrt(Ix.2+Iy.2);k=2;g=exp(-deltI./k);gx,gy=gradient(g)
6、;% 初始化圆,定义中心和半径center=floor(size(im)/2);radius = min(center)-8;u = init_u( imsize, center, radius);% 调用迭代函数filename = 3.bmp;m_name = filename( 1 : strfind( filename, . ) - 1 );num=400;u_new=die_dai(im,u,g,num,m_name);2、主要子程序(1)初始化水平集函数:function u = init_u( imsize, center, radius )% 初始化水平集m = imsize(
7、 1 ); n = imsize( 2 );u = zeros( imsize );for i = 1 : m;for j = 1 : n;distance = sqrt( sum( ( center - i, j ).2 ) );u( i, j ) = distance - radius;endend(2)迭代计算函数:function u=die_dai(im,u,g,num,m_name)m,n=size(u);gx,gy=gradient(g);u1=u;newpic=im; for i=2:m-1 for j=2:n-1if (u(i,j)*u(i+1,j)0u(i,j)=min(d
8、ist);elseu(i,j)=-min(dist);endendend五 实验结果分析以下为实验结果:Fig.1 Fig.2 Fig.3 Fig.4Fig.5 Fig.6 Fig.7 Fig.8Fig.9 Fig.10 Fig.11上面的图是每迭代 40 次选取一副图象,直至演化至稳态的结果(Fig.1 为调整大小并滤波后的图) 。注意:1曲线最终演化基本稳定在对象边界。2由于 c 选择过小,使曲线在对象凹陷部分(k0)继续运动的作用力()不足,但 c 过大又会使曲线穿透边界,因此 c 的选择要合适。guk3重新初始化间隔不要过大,会积累误差。4 中 k 值选取过大,导致函数 g 的下降速度
9、过于平缓,不能在对象边缘exp()r有效的停止。5 时间步长 delt_t 选取过大,就会影响效果。对不同的参数进行验证,实验结果表明各参数较为合理的值的范围:c:5 10delt_t: 0.05 0.1re_init:5 20k:0.0010.005六 实验小结:通过本次实验,可以看到 GAC(测地线活动轮廓模型)在图象分割中的理论运用。其数学思想是通过最小化一个能量泛函来实现的。这与前面所说的演化有些类似。只是在选取时采用什么样的函数来做不同。对于图象的分割,就是要找到边缘部分(即梯度变化在局部有最大值的地方) 。所以用 g(r)函数来做,使得在图像边缘处曲线停止演化,得到相应的水平集方法的 PDE,能够有效的进行分割。