智能计算导论课程设计

西安电子科技大学

智能计算导论课程实验报告SAR图像变化检测

SAR图像变化检测

1 引言

遥感变化检测是指通过对不同时期同一区域的遥感图像进行比较分析，根

据图像之间的差异得到我们所需要的地物或目标的变化信息。现代遥感技术的飞速发展为变化检测提供了一种便捷的途径，遥感数据成为变化检测的主要数据源。

与可见光和红外遥感相比，微波遥感具有无可比拟的优点：微波能穿透云雾、雨雪，具有全天候、全天时的工作能力。二，微波对地物有一定穿透能力。三，采用侧视方式成像，覆盖面积大。正是这些优点，使得SAR 图像日益成为变化

检测的重要数据源。

SAR 变化检测技术的需求日益广泛。目前，全球坏境变化加剧，城市急速发展，洪水、地震等自然灾害时有发生，这些都需要及时掌握相关动态信息，为相关决策部门提供支持，而SAR 的种种优点为快速响应提供了技术支持和应急保障。

2 定义

变化检测是指通过分析同一地区不同时间的图像，检测出该地区的地物随时间发生变化的信息

SAR图像的变化检测是指利用多时相获取的同一地表区域的SAR图像来确定和分析地表变化，能提供地物的空间展布及其变化的定性与定量信息

3 常用方法

本报告所用算法流程图

4 实验结果

5 程序

%initial clc

clear all close all

Ia=imread('2002.5.bmp');%read image

Ib=imread('2005.4.bmp');

Iag=Ia(:,:,1);%rgb2gray

Ibg=Ib(:,:,1);

%midfilt

Iam=medfilt2(Iag);

Ibm=medfilt2(Ibg);

%find difference

ia_double = double(Iam)+1; %uint8todouble

ia1 = ia_double/255; %unit

ib_double = double(Ibm)+1;

ib1 = ib_double/255;

di_image=di(ia1,ib1); %构造模糊差异

di_image1=uint8(di_image.*256);

di_image1=double(di_image1);

%FCM

bilateral_di_image=bilateral(di_image1,36,6); %双边滤波

bilateral_di_image1=uint8(bilateral_di_image);

bilateral_di_image=double(bilateral_di_image1);

fcm_image=fcm(bilateral_di_image); %聚类

[T]=mis(Re,fcm_image,N,L);

k=T;

disp(sprintf('¹²ÓÐ%d¸öÏñËØ£¬Îó²îÏñËØΪ%d',N*L,k));

p=k/(N*L)*100

disp(sprintf('Îó²îÂÊΪ%2.4f',p))

%display

figure(1);

subplot(2,3,1),imshow(Ia); title('原图a');

subplot(2,3,2),imshow(Iam); title('图a中值滤波图');

subplot(2,3,3),imshow(Ib); title('原图b');

subplot(2,3,4),imshow(Ibm); title('图b中值滤波图');

subplot(2,3,5),imshow(di_image); title('模糊差异图');

subplot(2,3,6),imshow(bilateral_di_image1); title('双边滤波图'); figure(2);

imshow(fcm_image); title('聚类图');

function [img1] = bilateral(img,winsize,sigma)

winsize = round( (winsize-1)/2 )*2 + 1;

[H,W] = size(img);

%Gaussian spatial filter

g_filter = fspecial('gaussian',winsize,winsize/6);

padnum = (winsize-1)/2;

A = padarray(img,[padnum padnum],'replicate','both');

img1 = zeros(size(img));

for jj = padnum+1:(padnum+1+H-1)

for kk = padnum+1:(padnum+1+W-1)

% Get a local neighborhood

imgwin = A(jj-padnum:jj+padnum,kk-padnum:kk+padnum);

% Find weights according to intensity diffs

Wwin = exp(-abs(imgwin - imgwin(padnum+1,padnum+1))/sigma^2);

% Find composite filter

newW = Wwin.*g_filter;

t = sum(sum(newW));

if(t>0)

newW = newW/t;

end

img1(jj-padnum,kk-padnum) = sum(sum(imgwin.*newW));

end

end

function di=di(IM1,IM2)

[sizeX,sizeY]=size(IM1);

A1=sizeX;

B1=sizeY;

[sizeX,sizeY]=size(IM2);

A2=sizeX;

B2=sizeY;

tmp1=zeros(A1,B1); %最小矩阵窗

tmp1=zeros(A1,B1); %最大矩阵窗

for i=1:1:A1

for j=1:1:B1

if IM1(i,j)>=IM2(i,j)

tmp1(i,j)=IM2(i,j);

tmp2(i,j)=IM1(i,j);

else tmp1(i,j)=IM1(i,j);

tmp2(i,j)=IM2(i,j);

end

end

end

mintmp=zeros(A1-2,B1-2);

maxtmp=zeros(A1-2,B1-2);

for i=1:1:A1