clear all
clc
lines=400;
samples=640;
% lines=1000;
% samples=1000;
N=2;
img=fopen('sample_BSQ','rb');
%img=fopen('11131','rb');
%for i=1:N
bi=fread(img,lines*samples,'uint8');
band_cov=reshape(bi,samples,lines);
band_cov2=band_cov';
band_uint8=uint8(band_cov2);
tif1=imadjust(band_uint8);
%mkdir('MATLAB_BSQ','tifbands')
%name=['MATLAB_BSQ\tifbands\tif',int2str(i),'.tif'];
%imwrite(tif1,'D:/Document/vc/matlab/图像分割/band1.tif','tif');
figure,imshow(tif1);
hold on;
%end
% bi=fread(img,lines*samples,'uint16');
% band_cov=reshape(bi,samples,lines);
% band_cov2=band_cov';
% band_uint8=uint16(band_cov2);
% tif2=imadjust(band_uint8);
%
% %mkdir('MATLAB_BSQ','tifbands')
% %name=['MATLAB_BSQ\tifbands\tif',int2str(i),'.tif'];
% %imwrite(tif2,'D:/Document/vc/matlab/图像分割/band2.tif','tif');
%
% figure,imshow(tif2);
% hold on;
fclose(img);
%% 提取边缘
%% 第一种方法:拉普拉斯算子
% g=zeros(lines-2,samples-2);
% for ii=2:lines-1
% for jj=2:samples-1
% g(ii-1,jj-1)=tif1(ii+1,jj)+tif1(ii-1,jj)+tif1(ii,jj+1)+tif1(ii,jj-1)-4*tif1(ii,jj);
% if g(ii-1,jj-1)>0
% g(ii-1,jj-1)=255;
% else
% g(ii-1,jj-1)=0;
% end
% end
% end
% figure,imshow(g);
% hold on;
%%
%% 第二种方法:马尔算子
% tif1r=reshape(tif1,lines*samples,1);
% tif1rd=double(tif1r);
% d=var(tif1rd);
% h=zeros(lines,samples);
% for ii=1:lines
% for jj=1:samples
% h(ii,jj)=exp(-(ii^2+jj^2)/(2*d));
% end
% end
% tif1d=double(tif1);
% g=h.*tif1d;
% g1=zeros(lines-2,samples-2);
% for ii=2:lines-1
% for jj=2:samples-1
% g1(ii-1,jj-1)=g(ii+1,jj)+g(ii-1,jj)+g(ii,jj+1)+g(ii,jj-1)-4*g(ii,jj);
% if g1(ii-1,jj-1)>0
% g1(ii-1,jj-1)=255;
% else
% g1(ii-1,jj-1)=0;
% end
% end
% end
% g2=uint8(g1);
% figure,imshow(g2);
%%
%% 第三种方法:Canny检测算子
k=1;
gx=zeros(lines,samples);
gy=zeros(lines,samples);
m=zeros(lines,samples);
a=zeros(lines,samples);
tif1d=double(tif1);
for ii=1:lines-1
for jj=1:samples-1
gx(ii,jj)=(tif1d(ii+1,jj)-tif1d(ii,jj)+tif1d(ii+1,jj+1)-tif1d(ii,jj+1))/2;
gy(ii,jj)=(tif1d(ii,jj)+tif1d(ii+1,jj)-tif1d(ii,jj+1)-tif1d(ii+1,jj+1))/2;
m(ii,jj)=sqrt(gx(ii,jj)^2+gy(ii,jj)^2);
a(ii,jj)=180*atan(gy(ii,jj)/gx(ii,jj))/pi;
end
end
%m1=uint8(m);
%figure,imshow(m1);
m1=zeros(lines,samples);
for ii=2:lines-1
for jj=2:samples-1
if (a(ii,jj)>67.5&&a(ii,jj)<90)||(a(ii,jj)>(-90)&&a(ii,jj)<(-67.5))
tmp1=m(ii,jj-1);
tmp2=m(ii,jj+1);
if m(ii,jj)>tmp1&&m(ii,jj)>tmp2
m1(ii,jj)=m(ii,jj);
else
m1(ii,jj)=0;
end
end
if a(ii,jj)>22.5&&a(ii,jj)<=67.5
tmp1=m(ii+1,jj-1);
tmp2=m(ii-1,jj+1);
if m(ii,jj)>tmp1&&m(ii,jj)>tmp2
m1(ii,jj)=m(ii,jj);
else
m1(ii,jj)=0;
end
end
if a(ii,jj)>(-22.5)&&a(ii,jj)<=22.5
tmp1=m(ii+1,jj);
tmp2=m(ii-1,jj);
if m(ii,jj)>tmp1&&m(ii,jj)>tmp2
m1(ii,jj)=m(ii,jj);
else
m1(ii,jj)=0;
end
end
if a(ii,jj)>=(-67.5)&&a(ii,jj)<=(-22.5)
tmp1=m(ii-1,jj-1);
tmp2=m(ii+1,jj+1);
if m(ii,jj)>tmp1&&m(ii,jj)>tmp2
m1(ii,jj)=m(ii,jj);
else
m1(ii,jj)=0;
end
end
end
end
m2=uint8(m1);
for ii=1:lines
for jj=1:samples
if m1(ii,jj)>(35)
m1(ii,jj)=m1(ii,jj);
else
m1(ii,jj)=0;
end
end
end
m2=uint8(m1);
figure,imshow(m2);
imwrite(m2,'D:/Document/vc/matlab/图像分割/bianyuan1.tif','tif');
%% 均值滤波
% tif2=zeros(lines,samples);
% for ii=1:lines-2
% for jj=1:samples-2
% tif2(ii,jj)=(tif1(ii,jj)+tif1(ii+1,jj)+tif1(ii+1,jj+1)+tif1(ii+1,jj+2)+tif1(ii+2,jj+1)+tif1(ii+2,jj+2)+tif1(ii,jj+2)+tif1(ii+2,jj)+tif1(ii,jj+1))/9;
% end
% end
% tif2u=uint8(tif2);
% figure,imshow(tif2u);
%% 中值滤波
% tif2=zeros(lines,samples);
% for i = 1:lines-2
% for j = 1:samples-2
% w= tif1(i:i+2,j:j+2);
% w1= sort(w);
% midum = w1(2,2,:);
% tif2(i,j) = midum;
% end
% end
% tif2u=uint8(tif2);
% figure,imshow(tif2u);