【图像分割】各向异性热扩散方程SAR图像分割【含Matlab源码 1879期】

news2024/11/23 3:53:57

⛄一、各向异性热扩散方程简介（具体理论见参考文献）

1 前言
图像分割是图像解译的一个关键步骥。图像分割属于图像解译的早期处理，目的是为更高水平的后续处理提供输人。因此图像分割的好坏对后续解译性能有较大的影响。相对于光学图像而言，图像的最大特点在于相干斑噪声的影响，它的存在使得图像表现为低信噪比，因此许多标准的图像分割算法应用于图像很难得到满意的效果。到目前为止，针对图像的特点，人们提出了许多图像分割算法上。这些方法可以归纳为两种分割思路首先对原始图像滤波，以降低斑点噪声的影响，然后采用与处理光学图像相似的方法进行分割，结合图像中象素的强度信息以及结构信息，研究在分割的过程中抑制斑点噪声影响的分割算法。相对于第二种思路，第一种思路带来的不利因素有对视数较大的图像分割可以取得较好的效果，但在噪声较强的情况下，为了达到好的去斑效果，会增大滤波的程度, 使得分割效果下降。而且降斑的客观全面评价标准不易确定 ,对分割产生的作用难以评估因此很难给出一个确定性的结论, 即到底对于不同质量的图像, 降斑到什么程度才能取得好的分割效果引人了计算复杂度 ,性能好的降斑算法时间可能很长,使得分割的速度很慢。因此 , 近年来关于图像分割的研究主要集中在第二种思路上,有代表性的方法如基于随机场的分割方法。这类方法分割精度很高,但其缺点在于分割的速度较慢,算法有可能陷人局部优化。同时值得注意的是文献中利用非线性热扩散方程进行图像的分割 ,该方法证明了非线性热扩散方程对于乘性噪声的抑制作用，但该方法需要已知图像的先验分布参数，本质上是一种有监督的分割方法 , 因而无法应用于实际的跟系统。本文采用第二种研究思路，提出了一种图像快速分割方法。该方法首先利用像素的强度信息，结合实际图像的统计分布，给出准则下最优分割，然后把该分割结果作为初始分割，采用各向异性热扩散方程对每类后验概率矩阵进行多尺度平滑，以滤除斑点噪声的影响

⛄二、部分源代码

function Flag = RSRADFilter(P_Class,p);
SIZE = size(P_Class);
SIZEX = SIZE(1,1);
SIZEY = SIZE(1,2);
landa = 8.8;rou = 2;
n = 12;q0 = 0.33333;%n = 20
P_Class0 = P_Class;
%初始化数组
I_gradR = zeros(SIZEX,SIZEY,2);I_gradL = zeros(SIZEX,SIZEY,2);
P_Class2 = zeros(SIZEX+2,SIZEY+2);
Cr2 = zeros(SIZEX+2,SIZEY+2);
% for i = 1:SIZEX
% for j = 1:SIZEY
% if P_Class(i,j) == 0;
% P_Class(i,j) = 1;
% end
% end
% end
%求图像I的R方向和L方向上的梯度
for time = 1:n
str = ‘已经完成’;
str = strcat(str,num2str(time/n100.0),‘%’);
disp(str);
q0_t(time) = q0exp(-1routime0.05);
% Delta_T = 5^{0.5*q0_t(time)}2;
%设定图像边界条件
P_Class2(2:SIZEX+1,2:SIZEY+1) = P_Class;
P_Class2(2:SIZEX+1,1) = P_Class(:,1);P_Class2(1,2:SIZEY+1) = P_Class(1,:);
P_Class2(2:SIZEX+1,SIZEY+2) = P_Class(:,SIZEY);P_Class2(SIZEX+2,2:SIZEY+1) = P_Class(SIZEX,:);
for i = 1:SIZEX
for j = 1:SIZEY
I_gradR(i,j,:) = [P_Class2(i+1+1,j+1)-P_Class2(i+1,j+1),P_Class2(i+1,j+1+1)-P_Class2(i+1,j+1)];
I_gradL(i,j,:) = [P_Class2(i+1,j+1)-P_Class2(i+1-1,j+1),P_Class2(i+1,j+1)-P_Class2(i+1,j-1+1)];
module_grad = I_gradR(i,j,1)^{2+I_gradR(i,j,2)}2+I_gradL(i,j,1)^{2+I_gradL(i,j,2)}2;
%计算Laplacian算子
I_Laplace(i,j) = P_Class2(i+1+1,j+1)+P_Class2(i+1,j+1+1)+P_Class2(i-1+1,j+1)+P_Class2(i+1,j-1+1)-4P_Class2(i+1,j+1);
%计算扩散系数
q(i,j) = ((1/4module_grad/(P_Class2(i+1,j+1)+0.0001)^{2-1/16*(I_Laplace(i,j)/(P_Class2(i+1,j+1)+0.0001))}2)/…
(1+1/4I_Laplace(i,j)/(P_Class2(i+1,j+1)+0.0001)+0.0001)²⁾0.5;
Cr(i,j) = p1/(2(1+(q(i,j)-q0_t(time))²⁾0.5);%2.53
% Cr(i,j) = 0.36exp(log2(1+(q(i,j)-q0_t(time))^{2))/((1+(q(i,j)-q0_t(time))}2)0.3010);
% Cr(i,j) = 2.361/((23.14)^{0.5*0.7)*exp(-1*(q(i,j)-q0_t(time)).}2/(20.7^2));% 1.88
% Cr(i,j) = 1.81/((23.14)^{0.66*0.4)*exp(-1*10*(q(i,j)-q0_t(time)).}2/(20.66^2));
% Cr(i,j) = 21/((23.14)^{1.66*0.4)*exp(-1*(q(i,j)-q0_t(time)).}2/(21.66^2));
end
end
%设定Cr的边界条件
Cr2(2:SIZEX+1,2:SIZEY+1) = Cr;
Cr2(2:SIZEX+1,1) = Cr(:,1);Cr2(1,2:SIZEY+1) = Cr(1,:);
Cr2(2:SIZEX+1,SIZEY+2) = Cr(:,SIZEY);Cr2(SIZEX+2,2:SIZEY+1) = Cr(SIZEX,:);
%计算散度
for i = 1:SIZEX
for j = 1:SIZEY
DIV(i,j) = Cr2(i+1+1,j+1)(P_Class2(i+1+1,j+1)-P_Class2(i+1,j+1))+Cr2(i+1,j+1)(P_Class2(i-1+1,j+1)-…
P_Class2(i+1,j+1))+Cr2(i+1,j+1+1)(P_Class2(i+1,j+1+1)-P_Class2(i+1,j+1))+Cr2(i+1,j+1)*(P_Class2(i+1,j-1+1)-P_Class2(i+1,j+1));
% end
% end