搜档网
当前位置:搜档网 › 遥感数字图像处理实习报告含Matlab处理代码

遥感数字图像处理实习报告含Matlab处理代码

遥感数字图像处理实习报告含Matlab处理代码
遥感数字图像处理实习报告含Matlab处理代码

辽宁工程技术大学

《数字图像处理》上机实习报告

教学单位辽宁工程技术大学

专业摄影测量与遥感

实习名称遥感数字图像处理

班级测绘研11-3班

学生姓名路聚峰

学号471120212

指导教师孙华生

实习1 读取BIP 、BIL、BSQ文件

一、实验目的

用Matlab读取BIP 、BIL、BSQ文件,并将结果显示出来。

遥感图像包括多个波段,有多种存储格式,但基本的通用格式有3种,即BSQ、BIL和BIP格式。通过这三种格式,遥感图像处理系统可以对不同传感器获取的图像数据进行转换。BSQ是像素按波段顺序依次排列的数据格式。BIL 格式中,像素先以行为单位块,在每个块内,按照波段顺序排列像素。BIP格式中,以像素为核心,像素的各个波段数据保存在一起,打破了像素空间位置的连续性。

用Matlab读取各个格式的遥感数据,是图像处理的前提条件,只有将图像读入Matlab工作空间,才能进行后续的图像处理工作。

二、算法描述

1.调用fopen函数用指定的方式打开文件。

2.在for循环中调用fread函数,用指定的格式读取各个像素。

3.用reshape函数,重置图像的行数列数。

4.用imadjust函数调整像素的范围,使其有一定对比度。

5.用imshow显示读取的图像。

三、Matlab源代码

1.读取BSQ的源代码:

clear all

clc

lines=400;

samples=640;

N=6;

img=fopen('D:\sample_BSQ','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);

tif=imadjust(band_uint8);

mkdir('D:\MATLAB','tifbands1')

name=['D:\MATLAB\tifbands1\tif',int2str(i),'.tif'];

imwrite(tif,name,'tif');

tilt=['波段',int2str(i)];

subplot(3,2,i),imshow(tif);title(tilt);

end

fclose(img);

2.读取BIP源代码

clear all

clc

lines=400;

samples=640;

N=6;

for i=1:N

img=fopen('D:\MATLAB\sample_BIP','rb');

b0=fread(img,i-1,'uint8');

b=fread(img,lines*samples,'uint8',(N-1));

band_cov=reshape(b,samples,lines);

band_cov2=band_cov';%×a??

band_uint8=uint8(band_cov2);

tif=imadjust(band_uint8);

mkdir('E:\MATLAB','tifbands')

name=['E:\MATLAB\tifbands\tif',int2str(i),'.tif'];

imwrite(tif,name,'tif'); %imwrite(A,filename,fmt)

tilt=['波段',int2str(i)];

subplot(3,2,i),imshow(tif);title(tilt);

fclose(img);

end

3.读取BIL的源代码

clear all

clc

lines=400;

samples=640;

N=6;

for i=1:N

bi=zeros(lines,samples);

for j=1:samples

img=fopen('D:\MATLAB\sample_BIL','rb');

bb=fread(img,(i-1)*640,'uint8');

b0=fread(img,1*(j-1),'uint8');

bandi_linej=fread(img,lines,'uint8',1*(N*samples-1));

fclose(img);

bi(:,j)=bandi_linej;

end

band_uint8=uint8(bi);

tif=imadjust(band_uint8);

mkdir('D:\MATLAB','tifbands')

name=['D:\MATLAB\tifbands\tif',int2str(i),'.tif'];

imwrite(tif,name,'tif');

tilt=['2¨??',int2str(i)];

subplot(3,2,i),imshow(tif);title(tilt);

end

四、运行结果

图1:读取文件的六个波段图

实习2 均值/中值滤波、边缘信息提取

一、实验目的与原理

各种图像滤波算子可以实现图像的增强,去噪,边缘提取等。图像增强的目的在于:1.采用一系列技术改善图像的视觉效果,提高图像的清晰度,2.将图像转换成一种更适合于人或机器进行分析处理的形式。它不是以图像保真度为原则,而是通过处理,设法有选择地突出便于人或机器分析某些感兴趣的信息,抑制一些无用的信息,以提高图像的使用价值。图像增强方法从增强的作用域出发,可分为空间域增强和频率域增强。空间域增强就是直接对图像像素灰度进行操作;频率域增强是对图像经傅里叶变换后的频谱成分进行操作,然后经傅里叶逆变换获得所需结果。

图像滤波可分为空间域滤波和频率域滤波,前者通过窗口或卷积核进行,它参照相邻像素改变单个像素的灰度值。后者对图像进行傅立叶变换,然后对频谱进行滤波。空间域图像滤波称为平滑和锐化,强调像素与其周围相邻像素的关系。去噪滤波为平滑滤波包括均值滤波和中值滤波。锐化滤波包括罗伯特梯度、索伯尔梯度、拉普拉斯算法、定向检测,用以提取线状地物和边缘。此实验用Matlab 采用各种滤波对图像进行了处理,处理结果如下:

二、算法描述

1.用imread读取图像文件,并用size获取图像的大小。

2.设计各种滤波算子。

3利用卷积公式对图像的没一个像素进行处理,得到滤波后的图像。

4.用imshow显示滤波后的图像。

三、Matlab源代码

1.均值滤波源码:

clear all

clc

img=imread('2.jpg');

[row,column,band]=size(img);

img0=double(img);

f11=1/9; f12=1/9; f13=1/9;

f21=1/9; f22=1/9; f23=1/9;

f31=1/9; f32=1/9; f33=1/9;

img1=[img0(:,1,:), img0(:,:,:), img0(:,column,:)];

img2=[img1(1,:,:); img1(:,:,:); img1(row,:,:)];

filtered=zeros(row,column,band);

for ii=1: row

for jj=1: column

filtered(ii,jj,:)=f11*img2(ii,jj,:) + f12*img2(ii,jj+1,:) + f13*img2(ii,jj+2,:)+ ...

f21*img2(ii+1,jj,:) + f22*img2(ii+1,jj+1,:) + f23*img2(ii+1,jj+2,:) + ...

f31*img2(ii+2,jj,:) + f32*img2(ii+2,jj+1,:) + f33*img2(ii+2,jj+2,:);

end

end

filtered1=uint8(filtered);

subplot(1,2,1),imshow(img);title('图1 原始RGB图像');

subplot(1,2,2),imshow(filtered1);title('图2 均值滤波后的图像');

imwrite(filtered1,'flower_filtered_mean.jpg');

2.边缘提取滤波源代码

clear all

img=imread('2.jpg');

[row,column,band]=size(img);

img0=double(img);

f11=1; f12=0; f13=-1;

f21=1; f22=0; f23=-1;

f31=1; f32=0; f33=-1;

img1=[img0(:,1,:), img0(:,:,:), img0(:,column,:)];

img2=[img1(1,:,:); img1(:,:,:); img1(row,:,:)];

filtered=zeros(row,column,band);

for ii=1: row

for jj=1: column

filtered(ii,jj,:)=f11*img2(ii,jj,:) + f12*img2(ii,jj+1,:) + f13*img2(ii,jj+2,:)+ ...

f21*img2(ii+1,jj,:) + f22*img2(ii+1,jj+1,:) + f23*img2(ii+1,jj+2,:) + ...

f31*img2(ii+2,jj,:) + f32*img2(ii+2,jj+1,:) + f33*img2(ii+2,jj+2,:);

end

end

filtered1=uint8(filtered);

subplot(1,2,1),imshow(img);title('图1 RGB原图像');

subplot(1,2,2),imshow(filtered1);title('图2 边缘提取后的图像');

imwrite(filtered1,'flower_filtered_edge.jpg');

四、运行结果

图1:原始RGB图像图2:均值滤波后的图像

图3:边缘提取后的图像

实习3 傅里叶变换、傅里叶逆变换,及频域滤波

一、实验目的

按照信号处理理论,根据滤除的频率特征,滤波有3种:1.低通滤波。低通滤波是对频率域的图像通过滤波器H(u,v)削弱或抑制高频部分而保留低频部分的滤波方法。由于图像上的噪声主要集中在高频部分,所以低通滤波可以起到压抑噪声的作用。同时,由于强调了低频成分,图像会变得比较平滑。2.高通滤波。高通滤波是对频率域的图像通过滤波器来突出图像的边缘和轮廓,进行图像锐化

的方法。3.带通滤波。仅保留指定频率范围的滤波,范围外的频率被阻止。

将空间域中的图像变换到频率域中进行计算。空间增强技术强调像元位置和像元之间的关系,但随着考虑的像元数目增多,计算的复杂度增加而且非常耗费计算运算时间,特别是当模板越来越大时,这种现象尤为明显。

频率域增强方法:

1.频率域平滑:保留图像的低频部分而抑制高频部分。

2.频率域锐化:保留图像的高频部分而削弱低频部分。

首先将空间域图像f(x,y)通过傅立叶变换为频率域图像F(u,v),然后选择合适的滤波器H(u,v)对F(u,v)的频谱成分进行增强得到图像G(u,v),再经过傅立叶逆变换将

G(u,v)变换到空间域,得到增强后的图像g(x,y)。

根据傅里叶变换的原理,用Matlab实现对图像的傅里叶变换,再设计各种频率滤波器,包括理想滤波器、巴特沃斯滤波器、指数滤波器等高通或低通滤波器,用这些滤波器对频率图像进行滤波,将得到的滤波后的频率图像经过傅里叶逆变换后得到想要的图像。本实验,对这些滤波器都进行了设计,处理结果如下图:

二、算法描述

1.读取RGB图像,并获得其某个波段。

2.调用fft2对某波段进行傅里叶变换,用fftshift频移函数,使得图像的低频部分移动到频谱的中心。

3.设置低通高通滤波器,对频谱图像进行滤波处理,去除其高频或低频部分。

4.用ifft2对频谱图像进行逆傅里叶变换,得到滤波后的图像。

5.对图像进行归一化,使像素值在0-255之间。

6. imshow显示频谱图像和滤波后的图像。

三、Matlab源代码

1.理想低通滤波源代码

clear all

a=imread('2.jpg');

[X,Y,Z]=size(a);

a1=a(:,:,1); b1= fft2(a1);

b2= fftshift(b1);

F=abs(b2);

h=zeros(X,Y);

cutoff=0.7; threshold=1-cutoff;

lowx=round(X/2-threshold*X/2);

upx=round(X/2+threshold*X/2);

lowy=round(Y/2-threshold*Y/2);

upy=round(Y/2+threshold*Y/2);

h(lowx:upx,lowy:upy)=1;

lowpass=b2.*h;

d0=ifft2(lowpass);

result=abs(d0);

result=uint8(result);

F1=log10(reshape(F, X*Y,1));

F2=uint8( ((F1 - min(F1))/(max(F1) - min(F1)))*255);

F3=reshape(F2, X, Y);

subplot(1,3,1), imshow(F3), title('Fig.1 ?傅里叶频谱');

subplot(1,3,2), imshow(h), title('Fig.2 理想低通滤波器');

subplot(1,3,3), imshow(result),title('Fig.3 ?理想低通滤波结果');

2.巴特沃斯低通滤波代码:

clear all

a=imread('2.jpg');

[X,Y,Z]=size(a);

a1=a(:,:,1);b1= fft2(a1);

b2= fftshift(b1);

F=abs(b2);

h=zeros(X,Y);

lowpass=zeros(X,Y);

n1=round(X/2);

n2=round(Y/2);

dmax=sqrt(n1^2 + n2^2);

cutoff=0.7;d0=(1-cutoff) * dmax;

n=3;

for x=1:X

for y=1:Y

d=sqrt((x-n1)^2+(y-n2)^2);

h(x,y)=1/(1+(d/d0)^(2*n));

lowpass(x,y)=b2(x,y).*h(x,y);

end

end

lp=ifft2(lowpass);

result=abs(lp);

result=uint8(result);

F1=log10(reshape(F, X*Y,1));

F2=uint8( ((F1 - min(F1))/(max(F1) - min(F1)))*255);

F3=reshape(F2, X, Y);

subplot(1,3,1), imshow(F3), title('Fig.1 ?傅里叶频谱');

subplot(1,3,2), imshow(h), title('Fig.2 巴特沃斯低通滤波器'); subplot(1,3,3), imshow(result),title('Fig.3 巴特沃斯低通滤波结果');

3.指数低通滤波的主要代码:

for x=1:X

for y=1:Y

d=sqrt((x-n1)^2+(y-n2)^2);

h(x,y)=1/(1+exp(-(d/d0))^(-n));

lowpass(x,y)=b2(x,y).*h(x,y);

end

end

4.指数高通滤波的主要代码:

for x=1:X

for y=1:Y

d=sqrt((x-n1)^2+(y-n2)^2);

h(x,y)=1/(1+exp(-(d/d0))^(n));

lowpass(x,y)=b2(x,y).*h(x,y);

end

end

四、运行结果

图1:傅里叶频谱

图2:理想低通滤波器图3:傅里叶低通滤波结果

图4:理想高通滤波器 图5:傅里叶高通滤波结果

图6:巴特沃斯低通滤波器 图7:巴特沃斯低通滤波结果

图8:巴特沃斯高通滤波器 图9:巴特沃斯高通滤波结果

图10:指数高通滤波器 图11:指数高通滤波结果

图12:指数低通滤波器 图13:指数低通滤波结果

实习4 主成分变换、主成分逆变换

一、实验目的与原理

遥感多光谱影像波段多,信息量大,对图像解译很有价值。但数据量太大,在图像处理计算时,也常常耗费大量的机时,占据大量的磁盘空间。实际上,一些波段的遥感数据之间都有不同程度的相关性,存在着数据冗余。多光谱变换方法可通过函数变换,达到保留主要信息、降低数据量、增强或提取有用信息的目的。其变换的本质是对遥感图像进行线性变换,使多光谱空间的坐标系按一定规律进行旋转。多光谱变换主要有两种变换:主成分变换和缨帽变换。

(1)主成分变换(K-L变换)的目的:

1.K-L变换是图像分析与模式识别中的重要工具,用于特征抽取,降低特

2.征数据的维数

3.K-L变换用于图像压缩时可以实现有损压缩和无损压缩

4.K-L变换后的N幅图像统计上互不相关,因此K-L变换可以去除图像数据

的相关性,提取主要信息。

用Matlab实现六个波段图像的主成分分析。求出六个波段的协方差阵,再求协方差阵的特征值、特征向量。只取特征值比较大的数,其占各个波段主要成分,特征值小的,大部分噪声等内容。

二、算法描述

主成分变换过程:

1.读取六个波段的图像。

2.将每个波段的所有像素存成一列,共六列,形成MN*6矩阵。

3.求出每个波段的平均像素值。

4.利用算出的平均像素值,求MN*6矩阵的协方差阵。

5.求协方差阵的特征值与特征向量。

6.将特征值按从大到小排列,特征向量与特征值对应。

主成分的逆变换过程:

(1)确定保留的波段数

(2)将去除的波段用0值替代

(3)将去噪后的PCA结果 SCORE*inv( COEFF)

(4)去中心化(各波段加上其均值)

(5)重新排列图像行列数

三、Matlab源代码

1.主成分分析代码

clear all

clc

num_band=6; %% 2¨??êy

b1=imread('tif1.tif');

b2=imread('tif2.tif');

b3=imread('tif3.tif');

b4=imread('tif4.tif');

b5=imread('tif5.tif');

b6=imread('tif6.tif');

[x,y]=size(b1);

img0=[reshape(b1,x*y,1),reshape(b2,x*y,1),

reshape(b3,x*y,1),reshape(b4,x*y,1),reshape(b5,x*y,1),reshape(b6,x*y, 1)];

img1=double(img0);

average=mean(img1);

img3=img1-repmat(average,x*y,1);

c=cov(img3);

[coff,lam] = eig(c);

lam1=lam(6:-1:1,6:-1:1);

coff2=coff(:,6:-1:1);

lam_val=diag(lam1);

lam_sum=sum(lam_val);

lam_n_sum=sum(lam_val(1:3));

prop=lam_n_sum/lam_sum;

result=img3*coff2;

result2=reshape(result,x,y,6);

result3=uint8(result2(:,:,1)-min(min(result2(:,:,1)))/max(max(result2 (:,:,1))-min(min(result2(:,:,1))))*255);

imshow(result3);

hold on

result3=uint8(result2(:,:,2)-min(min(result2(:,:,2)))/max(max(result2 (:,:,2))-min(min(result2(:,:,2))))*255);

figure,imshow(result3);

result3=uint8(result2(:,:,3)-min(min(result2(:,:,3)))/max(max(result2 (:,:,3))-min(min(result2(:,:,3))))*255);

figure,imshow(result3);

result3=uint8(result2(:,:,4)-min(min(result2(:,:,4)))/max(max(result2 (:,:,4))-min(min(result2(:,:,4))))*255);

figure,imshow(result3);

result3=uint8(result2(:,:,5)-min(min(result2(:,:,5)))/max(max(result2 (:,:,5))-min(min(result2(:,:,5))))*255);

figure,imshow(result3);

result3=uint8(result2(:,:,6)-min(min(result2(:,:,6)))/max(max(result2 (:,:,6))-min(min(result2(:,:,6))))*255);

figure,imshow(result3);

result2(:,:,4:6)=0;

result_new=reshape(result2,x*y,6);

result_ipca=result_new*inv(coff2);

result_ipcaz2=result_ipca+repmat(average,x*y,1);

result_ipcaz3=reshape(result_ipcaz2,x,y,6);

figure,imshow(uint8(result_ipcaz3(:,:,1)));

C = cat(3, uint8(result_ipcaz3(:,:,3)),

uint8(result_ipcaz3(:,:,4)),uint8(result_ipcaz3(:,:,5)));

figure,imshow(C);

四、运行结果

图1:第一主分量图2:第二主分量

图3 :第三主分量图4:第四主分量

图5:第五主分量图6:第六主分量

图7:逆变换后的图像图8:逆变换后的彩色合成图像由图可以看出,图像的主要成分主要集中在前三个主分量图像上,后三个图像大部分都是噪声。通过逆变换恢复原来的图像,这样,经过K-L变换后的前三个主分量,信噪比大,噪声相对小,突出了主要信息,达到了图像增强的目的。另外,经过K-L变换后,第一或前二或前三个主分量已经包含了绝大多数的地物信息,足够分析使用,同时数据量大大减少。应用中,常常只取前三个主分量作假彩色合成(如图8),数据量可减少到43%,既实现了数据压缩,也可作为分类前的特征选择。

实习5 缨帽变换

一、实验目的

缨帽变换(K-T变换)的目的:该变换也是一种坐标空间发生旋转的线性组合变换,但旋转后的坐标轴不是指向主成分方向,而是指向与地物特别是和植被生成以及土壤有密切关系的方向

二、算法描述

1.读取六个波段的图像。

2.将每个波段的所有像素排成一列,共六列组成一个矩阵。

3.设计转换矩阵。

4.图像矩阵与转换矩阵相乘,得到六个相互垂直的分量。

5.逐一显示各个分量。

三、Matlab源代码

clear all

clc

b1=imread('tif1.tif');

b2=imread('tif2.tif');

b3=imread('tif3.tif');

b4=imread('tif4.tif');

b5=imread('tif5.tif');

b6=imread('tif6.tif');

[x,y]=size(b1);

img0=[reshape(b1,x*y,1),reshape(b2,x*y,1),

reshape(b3,x*y,1),reshape(b4,x*y,1),reshape(b5,x*y,1),reshape(b6,x*y, 1)];

img0=double(img0);

KT_tm = [0.3561, 0.3972, 0.3904, 0.6966, 0.2286, 0.1596;

-0.3344,-0.3544,-0.4556, 0.6966,-0.0242,-0.2630;

0.2626, 0.2141, 0.0926, 0.0656,-0.7629,-0.5388;

0.0805,-0.0498, 0.1950,-0.1327, 0.5752,-0.7775;

-0.7252,-0.0202, 0.6683, 0.0631,-0.1494,-0.0274;

0.4000,-0.8172, 0.3832, 0.0602,-0.1095, 0.0985];

img1 = img0 * KT_tm';

min0 = min(img1);

max0 = max(img1);

dm = max0 - min0;

img11 = img1 - repmat(min0,x*y,1);

img2 = (img11./repmat(dm,x*y,1))*255;

image = reshape(img2,x,y,6);

figure,imshow(uint8(image(:,:,1))),title('KT变换第一分量:亮度'); figure,imshow(uint8(image(:,:,2))),title('KT变换的第二分量:绿度'); figure,imshow(uint8(image(:,:,3))),title('KT变换的第三分量:湿度'); figure

imshow(uint8(image(:,:,4))),title('KT变换的第四分量:黄度指数及噪声');四、运行结果

图9:第一分量:亮度图 图10:第二分量:绿度

图11:第三分量:湿度 图12:第四分量:黄度指数及噪声 用Matlab 编写代码实现缨帽变换,其基本思想就是:采用变换矩阵,使图像指向与地物特别是和植被生长以及土壤有密切关系的方向。变换后

123456[,,,,,]Y y y y y y y ,这六个分量相互垂直,前三个分量具有明确的地物意义:y1为亮度,实际上是TM 六个波段的加权和,反映出图像总体的反射值;y2为绿度,y3为湿度。

实习6 图像分类

(最小距离分类、最大似然分类、K 均值分类)

一、实验目的与原理

图像分类的目的是将图像中每个像元根据其不同波段的光谱亮度、空间结构特征或者其他信息,按照某种规则或算法划分为不同的类别。遥感图像分类就是对地球表面及其环境在遥感图像上的信息进行识别和分类,从而达到识别图像信息所对应的实际地物,提取所需地物信息的目的。

遥感图像分类的原理:由于地物的成分、性质、分布情况的复杂程度和成像条件不同,以及一个像元或瞬时试场里往往有两种或多种地物的情况,即混合像元,使得同类地物的特征向量也不尽相同,而且使得不同地物类型的特征向量之间的差别也不都是显著的。遥感图像的光谱特征通常是以地物在多光谱图像上的灰度体现出来的,即不同地物在同一波段图像上表现的灰度一般互不相同;同时,不同地物在多个波段上的灰度呈现规律也不相同,这就构成了我们在图像上区分不同地物的物理依据。

遥感图像传统的计算机分类方式有监督分类和非监督分类两种。监督分类包括最小距离分类法、平行多面体分类法、最大似然分类法。非监督分类包括K-

均值聚类法、ISODATA分类法。本实验将用Matlab实现监督分类的最小距离分类和最大似然分类,以及非监督分类中的K-均值分类。

三、算法描述

(1)最小距离分类算法描述

1.读取各个波段文件,并加载各类地物在各个波段的均值文件。

2.计算每个像元向量到各类地物的距离。

3.比较像元到各类地物的距离,取最小的距离,并将像元归为该类。

4.为各类赋予一定的颜色。

(2)最大似然分类算法描述:

1.读取BIP文件,将每个波段读为一列。

2.加载各类地物在各个波段的均值和方差

3.对每个像元向量计算属于各个类的概率,比较概率的大小,概率最大的将像元归为该类。

4.对每个类赋予一定的颜色。

5.显示分类后的图像。

(3)K-means分类算法描述:

1.为每个聚类确定一个初始的聚类中心(一定要在最小值和最大值之间,这一步非常重要)

2.将样本集中的每一个样本按照最小距离原则,分配到k个聚类中的某一个3.使用每个聚类中所有样本的均值作为新的聚类中心

4.如果聚类中心有变化则重复2、3步直到聚类中心不再变化为止

5.最后得到的k个聚类中心就是聚类的结果

四、Matlab源代码

(1)最小距离分类代码:

clear all

clc

x = 400;

y = 640;

N = 6;

image = zeros(x*y,N);

for i = 1:N

img=fopen('sample_BIP','rb');

b0 = fread(img,i-1,'uint8');

image(:,i) = fread(img,x*y,'uint8',N-1);

fclose(img);

end

load crop.txt;

load forest.txt;

load water.txt;

load soil1.txt;

load soil2.txt;

load soil3.txt;

dis = zeros(x*y,N);

dis(:,1) = sqrt(sum((image-repmat(crop(1,:),x*y,1)).^2,2));

dis(:,2) = sqrt(sum((image-repmat(forest(1,:),x*y,1)).^2,2)); dis(:,3) = sqrt(sum((image-repmat(water(1,:),x*y,1)).^2,2)); dis(:,4) = sqrt(sum((image-repmat(soil1(1,:),x*y,1)).^2,2)); dis(:,5) = sqrt(sum((image-repmat(soil2(1,:),x*y,1)).^2,2)); dis(:,6) = sqrt(sum((image-repmat(soil3(1,:),x*y,1)).^2,2)); min_dis = min(dis,[],2);

classes = zeros(x*y,3);

for i = 1:x*y

ind = find(dis(i,:)==min_dis(i));

switch ind

case 1

classes(i,:) = [255,0,0];

case 2

classes(i,:) = [0,255,0];

case 3

classes(i,:) = [0,0,255];

case 4

classes(i,:) = [255,255,0];

case 5

classes(i,:) = [255,0,255];

case 6

classes(i,:) = [0,255,255];

end

end

img2 = reshape(classes,y,x,3);

RGB_classes = cat(3,img2(:,:,1)',img2(:,:,2)',img2(:,:,3)'); imshow(uint8(RGB_classes));

title('最小距离分类');

(3)最大似然分类代码:

clear all

clc

x = 400;

y = 640;

N = 6;

image = zeros(x*y,N);

for i = 1:N

img=fopen('sample_BIP','rb');

b0 = fread(img,i-1,'uint8');

image(:,i) = fread(img,x*y,'uint8',N-1);

fclose(img);

end

load crop.txt;

load forest.txt;

load water.txt;

load soil1.txt;

load soil2.txt;

load soil3.txt;

all=cat(3,crop,forest,water,soil1,soil2,soil3);

cov = zeros(6,6,6);

mean = zeros(1,6,6);

for i = 1:6

mean(1,:,i) = all(1,:,i);

cov(:,:,i)=diag(all(2,:,i));

end

dk = zeros(x*y,6);

for i = 1:x*y

for j = 1:6

dk(i,j) = -log(det(cov(:,:,j)))*0.5 - 0.5*(image(i,:) - mean(1,:,j))/cov(:,:,j) * (image(i,:) - mean(1,:,j))';

end

end

dk_max = max(dk,[],2);

classes = zeros(x*y,3);

for i = 1:x*y

ind = find(dk(i,:)==dk_max(i));

switch ind

case 1

classes(i,:) = [255,0,0];

case 2

classes(i,:) = [0,255,0];

case 3

classes(i,:) = [0,0,255];

case 4

classes(i,:) = [255,255,0];

case 5

classes(i,:) = [255,0,255];

case 6

classes(i,:) = [0,255,255];

end

end

img2 = reshape(classes,y,x,3);

RGB_classes = cat(3,img2(:,:,1)',img2(:,:,2)',img2(:,:,3)'); imshow(uint8(RGB_classes));

title('最大似然分类');

(3)K-means分类代码:

clear all

clc

x = 400;

y = 640;

N = 6;

image = zeros(x*y,N);

for i = 1:N

img=fopen('sample_BIP','rb');

b0 = fread(img,i-1,'uint8');

image(:,i) = fread(img,x*y,'uint8',N-1);

fclose(img);

end

load ini_center.txt;

center = zeros(size(ini_center));

threshold = 1e-3;

max_d = 1;

while max_d > threshold

dis = zeros(x*y,N);

dis(:,1) = sqrt(sum((image-repmat(ini_center(1,:),x*y,1)).^2,2));

dis(:,2) = sqrt(sum((image-repmat(ini_center(2,:),x*y,1)).^2,2)); dis(:,3) = sqrt(sum((image-repmat(ini_center(3,:),x*y,1)).^2,2));

dis(:,4) = sqrt(sum((image-repmat(ini_center(4,:),x*y,1)).^2,2));

dis(:,5) = sqrt(sum((image-repmat(ini_center(5,:),x*y,1)).^2,2));

dis(:,6) = sqrt(sum((image-repmat(ini_center(6,:),x*y,1)).^2,2));

min_dis = min(dis,[],2);

classes = zeros(x*y,3);

crop = zeros(1,6);

forest = zeros(1,6);

water = zeros(1,6);

soil1 = zeros(1,6);

soil2 = zeros(1,6);

soil3 = zeros(1,6);

n_crop = 0;

n_forest = 0;

n_water = 0;

n_soil1 = 0;

n_soil2 = 0;

n_soil3 = 0;

for i = 1:x*y

ind = find(dis(i,:) == min_dis(i));

switch ind(1,1)

case 1

crop = crop + image(i,:);

n_crop = n_crop + 1;

classes(i,:) = [255,0,0];

case 2

forest = forest + image(i,:);

n_forest = n_forest + 1;

classes(i,:) = [0,255,0];

case 3

water = water + image(i,:);

n_water = n_water + 1;

classes(i,:) = [0,0,255];

case 4

soil1 = soil1 + image(i,:);

n_soil1 = n_soil1 + 1;

classes(i,:) = [255,255,0];

case 5

soil2 = soil2 + image(i,:);

n_soil2 = n_soil2 + 1;

classes(i,:) = [255,0,255];

case 6

soil3 = soil3 + image(i,:);

n_soil3 = n_soil3 + 1;

classes(i,:) = [0,255,255];

end

end

center(1,:) = crop./n_crop;

center(2,:) = forest./n_forest;

center(3,:) = water./n_water;

center(4,:) = soil1./n_soil1;

center(5,:) = soil2./n_soil2;

center(6,:) = soil3./n_soil3;

max_d = max(max(abs(center-ini_center)));

ini_center = center;

end

img2 = reshape(classes,y,x,3);

RGB_classes = cat(3,img2(:,:,1)',img2(:,:,2)',img2(:,:,3)'); imshow(uint8(RGB_classes)),title('kmeans·?àà?á1?');

四、运行结果

图1:最小距离分类

遥感课程理解与应用读书报告

遥感图像理解与应用读书报告 一、概述 遥感主要是为地学研究提供数据源,但由于遥感数据获取的是地物电磁辐射信息,不能直接用于各类地学研究,需要先对遥感数据进行处理。对遥感数据的处理,包括遥感数据预处理、遥感数据理解。预处理包括对数据几何校正和图像配准、图像融合、图像镶嵌与裁剪、大气校正,使遥感数据与地理上的空间点相对应,包括位置对应以及位置辐射值的对应。而对遥感图像的理解则是根据遥感数据中的辐射值确定不同位置的目标类型、属性等信息,使数据符合人的空间认知。 遥感图像理解这门课程首先介绍遥感原理,它是通过获取和分析目标的电磁辐射信息,来了解目标的其他属性信息的科学和技术。然后介绍遥感图像,包括获取手段(直接获取,间接由现有图件或照片获取)、显示方式(真彩色、假彩色、伪彩色)、图像属性(四种分辨率)等。再对遥感图像中的目标做介绍,指出遥感最重要目的就是获取遥感图像目标及其语义解释。遥感图像目标表现于图像中的特征包括图像可视特征、图像参数特征和光谱特征,信息丰富目标识别能力增强。其中图像可视特征包括目标的亮度、颜色、纹理、边沿、轮廓、形态、大小等;图像参数特征是经过计算后得到的用于描述图像特征的各种参数,如图像灰度的均值、方差,图像的比值,图像的协方差、各阶矩,图像在变换域中的频谱等;图像光谱特征由各个波段的光谱值决定,包括平均光谱值大小、光谱曲线的变化趋势和光谱曲线中对地物信息具有标示性意义的一些几何参数,如波峰、波谷、斜率等。在对遥感图像理解中,主要针对这些信息确定图像目标类型、属性等信息。 二、对遥感图像的理解 第四章对遥感图像的理解是重点内容。针对地学应用的遥感图像中的目标是各种地理客体,因而这里的遥感图像理解也就是遥感图像的地学理解。地学遥感图像理解则除了包括目标的几何关系、目标类别外,更重要的是理解目标的性质、性状、数量特征等。其内涵主要有:目标类别、地物空间关系、目标的性状(物理、化学、生物参数)、目标的数量特征。 遥感图像理解包括图像处理、图像分析和图像理解三个层次的内容。图像处理包括图像纠正(也叫数据预处理)和图像增强,纠正图像的几何误差和辐射误差,并突出所感兴趣的信息。图像分析包括边缘检测、图像分类、空间分析,提取感兴趣的目标和信息,对图形空间信息进行综合分析。课程主讲的图像理解侧重于计算机解译,从图像特征或光谱特征出发,有很多不同的方法。 从图像特征出发,图像特征理解步骤为图像预处理、目标检测、目标解释,主要是以图像中的目标为理解单元,一般多用于高分辨率遥感图像中。从光谱特征出发,侧重光谱特征的理解包括数据预处理、波段选择、图像分类三步,一般用于低分辨率遥感图像中。对遥感数据具体采用什么理解方法视具体情况而定,理想的方式是将图像特征与光谱特征结合进行。 三、侧重图像特征的图像理解 图像特征包括图像色调、颜色、阴影、形状、纹理、大小、位置、图型、相关布局、参数特征等。基于图像目标检测的策略,用图像区域表达图像中的目标,其中点目标和线目标

遥感实验报告

重庆交通大学 学生实验报告 实验课程名称遥感原理与应用 开课实验室测量与空间信息处理实验室 学院 2013 年级测绘工程专业 1班学生姓名刘文洋 学号 631301040126 开课时间 2015 至 2016 学年第 1 学期

目录 实验一 ENVI 视窗的基本操作 (2) 实验二遥感图像的几何校正 (4) 实验三遥感图像的增强处理 (8) 实验四遥感图像的变换 (12) 实验五遥感信息的融合 (15) 实验六遥感图像分类 --- 监督分类 (17) 实验七遥感图像分类 --- 非监督分类 (19) 实验八遥感图像分类后处理 (22)

实验一ENVI 视窗的基本操作 一、实验目的 初步了解目前主流的遥感图象处理软件 ENVI 的主要功能模块,在此基础上,掌握视窗操作模块的功能和操作技能,为遥感图像的几何校正等后续实习奠定基础。 二、实验内容 视窗功能介绍;文件菜单操作;显示数据;裁剪数据;合并波段 三、实验步骤 1、首先打开ENVI4.7软件,看见的只有菜单栏,如图所示: 2、打开每个下拉菜单浏览其下拉栏中都有哪些功能,比如:我们如果需要打开遥感文件,则可以选择File下的打开功能open image file,打开遥感图像如下图:

裁剪数据打开basic tools的resize data功能,如果需要对图像进行一系列处理,可以利用Transform,Classification等功能进行操作,在后续实验中我们也会用到其中的一些功能进行图像的一系列操作,到时候在详细叙述。 3、再熟悉了ENVI4.7的一些基本知识后我们可以简单地操作下,比如对一组数据分别用Gray Scale和Load RGB导入,看看两幅图的区别以及各自的优缺点。 四、实验结果分析 在这次的实验中,我们简单的熟悉了下ENVI4.7的一些功能,发现它是可以对遥感图像进行图像几何纠正,直方图均衡,监督分类,非监督分类等一系列操作,为我们后续利用软件对遥感图像处理打下了基础。

数字图像处理_旋转与幅度谱(含MATLAB代码)

数字图像处理实验一 15生医 一、实验内容 产生右图所示图像 f1(m,n),其中图像大小为256 ×256,中间亮条为128×32,暗处=0,亮处=100。 对其进行FFT: ①同屏显示原图f1(m,n)和FFT(f1)的幅度谱图; ②若令f2(m,n)=(-1)^(m+n)f1(m,n),重复 以上过程,比较二者幅度谱的异同,简述理由; ③若将f2(m,n)顺时针旋转90度得到f3(m,n),试显示FFT(f3)的 幅度谱,并与FFT(f2)的幅度谱进行比较; ④若将f1(m,n) 顺时针旋转90度得到f4(m,n),令f5(m,n) = f1(m,n) + f4(m,n),试显示FFT(f5)的幅度谱,指出其与 FFT(f1)和FFT(f4)的关系; ⑤若令f6(m,n)=f2(m,n)+f3(m,n),试显示FFT(f6)的幅度谱,并指出其与 FFT(f2)和FFT(f3)的关系,比较FFT(f6)和FFT(f5)的幅度谱。 二、运行环境 MATLAB R2014a 三、运行结果及分析 1.同屏显示原图f1(m,n)和FFT(f1)的幅度谱图:

50100150200250 100150200250 50100150200250 100150200250 2.令f2(m,n)=(-1)^(m+n )f1(m,n),对其进行FFT ,比较f2与f1幅度谱的异同,简述理由: 50100150200250 100150200250 50100150200250 100150200250 异同及理由:①空域:f2由于前边乘了系数(-1)^(m+n ),导致灰度值有正有负,而在MATLAB 的imshow 函数中默认把负值变为0(有些情况是取反),所以形成了如左图所示的黑白花纹。②频域:FFT(2)

遥感数字图像处理实习1

(1)以多波段组合方式将GeoTIFF格式的白银市TM原始数据转换为ENVI Standard 格式: 利用Basic Tools/Layer Stacking弹出对话框然后Import File,弹出对话框,导入GeoTIFF格式的TM原始数据,选择波段1、2、3、4、5和7, 点击OK,利用Choose选择输出路径及文件名,同时可以利用Reorder Files对输入的文件根据自己的需要进行调换顺序,点击OK输出ENVI Standard格式的数据。 (2)查询并记录影像文件的基本信息、投影信息,以及各个波段直方图信息,然后编辑头文件: 利用Basic Tools/Resize Data弹出对话框里面选择要查看的影像,左 边会出现其基本信息,如图所示:也有投影信息,既可以用来看单波段的也可以看合成后整个影像的信息。在对话框下,合成影像的名字上右击,选择Quick Statistics弹出对话框,在此对话框中点击Select Plo下拉菜单,选择单波段或者多波段的直方图,相应的对话框中会出现直方图(在结果与分析中记录),还可以右击选择edit修改横、纵坐标的单位。 同样的在合成影像的名字上右击,选择Edit Head,弹出对话框

然后点击Edit Attributes/Band Name弹出对话框,选中波段输入修改 后名字,点击OK即可进行波段名字的修改。点击Edit Attributes/Wavelengths弹出进行相应的波长的修改。 (3)在View视窗中,利用影像缩小、放大、漫游工具识别影像中的土地利用/土地覆盖类型: 可以结合当地的google earth上高分辨率的遥感影像,进行识别,利用Viewer视窗下Tools/SPEAR/Google Earth/Jump to Location可以在google earth上显示View主视窗中相应选中地物对应的位置。 (4)利用Viewer视窗打开影像,分别选取4、3、2和7、4、2波段组合进行假彩色合成,观察实习内容中所要求地物的色调变化: 利用File/Open Image File,选择第1步合成的ENVI Standard 格式的数据,弹出对话框,在其中选择RGB Color,将R、G、B分别设为4、3、2波段,点击Load Band,在Viewer#1中出现了4、3、2波段组合的假彩色图像,再在此窗口中,点击Display/New Display,弹出Viewer#2,选择RGB Color,将R、G、B分别设为7、4、2波段,点击Load Band,在Viewer#2中出现了7、4、2波段组合的假彩色图像,在Viewer窗口中右击选择Link Displays,弹出对话框,点击OK,可以把两个窗口中同一位置进行连接起来, 即其中一个窗口放大、缩小、漫游到某个位置,另外一个也跟着漫游到其相对应的位置。这样可以进行地物色调变化的对比。 (5)提取6种地物在不同波段的数值(Digital Number,DN),做光谱剖面图: 在Viewer视窗中Tools/Profile/Z Profile(Spectrum)弹出对话框,在其 Options下拉菜单中勾选Plot Key,对话框中出现了Viewer视窗中选中的目标地物的X,Y坐标,然后勾选Collect Spectra,鼠标箭头变为十字箭头,在目标地物中取九个点(本来图上就有一个,总共是十个点),然后在选择File/Save Plot As/ASCII弹出对话框 ,点击Select All Items,利用Choose选择输出路径和文件名,点击 OK,将其保存为.txt格式。选六种地物,重复以上操作,提取不同波段的数值(Digital Number,DN)。将.txt格式的文件用excel打开,然后用插入函数中的average函数求出每种地物的平均DN值,然后做出光谱剖面(光谱图如结果与分析中所示)。 (6)使用Excel制作6种地物的样本特征光谱统计表: 在Excel中分别使用插入函数中的AVERAGE、VAR、STDEV、MAX和MIN函数求出各地物样本DN值在各个波段的平均值、方差、标准差、最大值和最小值。然后,在07版Excel 的“Microsoft Office 按钮”,单击“Excel 选项”。“加载项”,然后在“管理”框中,选择“Excel 加载项”,单击“转到”弹出“加载宏”,在弹出来的对话框中选择“分析工具库”,并点击确定。然后从“工具”中找到“数据分析”,从“数据分析”对话框中选择“协方差”,并导入某种地物需求协方差的数据区域并选择“逐行”进行,最后选择数据输出区域并确定,则可得该地物的协方差矩阵。同理,在从“数据分析”对话框中选择“相关系数”,进行相应操作,可求得相关系数矩阵。(在结果与分析中附有个地物的样本特征光谱统计表)(7)制作散点图: 在Excel中,打开6种地物的样本DN数据(5步骤产生的),选择band2和band4做散

中国地质大学遥感图像处理上机实习报告

遥感图像处理课程实习报告 学生姓名:王蜀越 班学号: 学号: 指导教师:王红平、许凯 中国地质大学信息工程学院 2017年7月1日

目录 目录 ............................................................................................................................................... - 1 - 实习一:影像融合........................................................................................................................ - 2 - 1.1【实习目的】 (2) 1.2【实习步骤】 (2) 1.3【实习过程】 (2) 实习二:几何校正........................................................................................................................ - 6 - 2.1【实习内容】 (6) 2.2【实习步骤】 (6) 2.3【实习过程】 (6) 实习三:影像分类(一).......................................................................................................... - 10 - 3.1【实习内容】 (10) 3.2【实习步骤】 (10) 3.3【实习过程】 (10) 实习四:影像分类(二).......................................................................................................... - 14 - 4.1【实习内容】 (14) 4.2【实习步骤】 (14) 4.3【实习过程】 (14) 心得与感想 ................................................................................................................................. - 18 -

遥感图像处理实验

哈尔滨工业大学 遥感图像处理及遥感系统仿真 实验报告 项目名称:《遥感图像处理及遥感系统仿真创新》 姓名:蒋国韬 学号:24 院系:电子与信息工程学院 专业:遥感科学与技术 指导教师:胡悦 时间:2017年7月

实验一:遥感数字图像的增强 一、实验目的: 利用一幅城市多光谱遥感图像,分析其直方图,并利用对比度增强和去相关拉伸方法对遥感图像进行增强。 二、实验过程: 1.用multibandread语句读取一幅多光谱遥感图像(7波段,512x512图像)的可 见1,2,3波段(分别对应R,G,B层); 2.显示真彩色图像; 3.通过研究直方图(imhist),分析直接显示的真彩色图像效果差的原因;

4.利用对比度增强方法对真彩色图像进行增强(imadjust,stretchlim); 5.画出对比度增强后的图像红色波段的直方图;

6.利用Decorrelation去相关拉伸方法(decorrstretch)对图像进行增强;

7.显示两种图像增强方法的结果图像。

三、实验分析: (1)高光谱影像由于含有近百个波段,用matlab自带的图像读写函数imread和imwrite往往不能直接操作,利用matlab函数库中的multibandred函数,可以读取多波段二进制图像。512×512为像素点,7位波段数,bil为图像数组的保存格式,uint8=>uint8为转换到matlab 的格式,[3 2 1]的波段分别对应RGB三种颜色。 (2)直接观察真彩复合图像发现,图像的对比度非常低,色彩不均匀。通过观察红绿蓝三色的波段直方图,可以观察到数据集中到很小的一段可用动态范围内,这是真彩色复合图像显得阴暗的原因之一。另外,根据三种颜色的三维散点图,如下

matlab图像处理代码

附录 MATLAB图像处理命令  1.applylut  功能: 在二进制图像中利用lookup表进行边沿操作。 语法: A = applylut(BW,lut) 举例 lut = makelut('sum(x(:)) == 4',2); BW1 = imread('text.tif'); BW2 = applylut(BW1,lut); imshow(BW1) figure, imshow(BW2) 相关命令: makelut 2.bestblk  功能: 确定进行块操作的块大小。 语法: siz = bestblk([m n],k) [mb,nb] = bestblk([m n],k) 举例 siz = bestblk([640 800],72) siz = 64 50 相关命令: blkproc 3.blkproc  功能:

MATLAB 高级应用——图形及影像处理 320 实现图像的显式块操作。 语法: B = blkproc(A,[m n],fun) B = blkproc(A,[m n],fun,P1,P2,...) B = blkproc(A,[m n],[mborder nborder],fun,...) B = blkproc(A,'indexed',...) 举例 I = imread('alumgrns.tif'); I2 = blkproc(I,[8 8],'std2(x)*ones(size(x))'); imshow(I) figure, imshow(I2,[]); 相关命令: colfilt, nlfilter,inline 4.brighten  功能: 增加或降低颜色映像表的亮度。 语法: brighten(beta) newmap = brighten(beta) newmap = brighten(map,beta) brighten(fig,beta) 相关命令: imadjust, rgbplot 5.bwarea  功能: 计算二进制图像对象的面积。 语法: total = bwarea(BW) 举例 BW = imread('circles.tif'); imshow(BW);

遥感读书报告

东南大学交通学院测绘工程系 遥感读书报告 专业:测绘工程 班级:213 学号:21311112 姓名:白金睿 老师:戚浩平 日期:2013年12月

首先纵观遥感这本书,我们可以先粗略的先把它分为35个小的标题,然后总结下来,之后再细分其中的重点,之后再详细说一说他们中的难点,就其中的一些我比较感兴趣的公式做一些推导解释,在上理论课与实验课的时候我慢慢发现RS的强大与其特有的魅力,ERDAS,ARCGIS,ARCMAP,这些软件即强大又有趣,让我来带领大家,纵游遥感的海洋吧~ 1.遥感技术系统的组成 被测目标的信息特征、信息的火枪、信息的传输与记录、信息的处理和信息的应用。2.遥感的类型 1)按遥感平台分为地面遥感、航空遥感、航天遥感; 2)按工作方式分为主动遥感和被动遥感; 3)按探测波段分为:紫外遥感(0.3-0.4);可见光(0.4-0.7);红外(0.7-14mm); 微波(0.1-100cm)等。 3.遥感技术的特点 大面积的同步观测、时效性、数据的综合性和可比性、经济性、局限性。 4.电磁波的主要参数 1)波长(Wavelength):指波在一个振动周期内传播的距离。即沿波的传播方向,两个相邻的同相位点(如波峰或波谷)间的距离。 2)周期:波前进一个波长那样距离所需的时间。 3)频率(frequency):指单位时间内,完成振动或振荡的次数或周期(T),用V示。 注:一般可用波长或频率来描述或定义电磁波谱的范围。在可见光——红外遥感中多用波长,在微波遥感中多用频率。 4)振幅(Amplitude):表示电场振动的强度。它被定义为振动物理量偏离平衡位置的最大位移,即每个波峰的高度。 5)电磁波谱:将各种电磁波在真空中的波长按其长短,依次排列制成的图表。5.常用电磁波波段特性 1)紫外线(UV):0.01-0.4μm,碳酸盐岩分布、水面油污染; 2)可见光:0.4-0.76 μm,鉴别物质特征的主要波段;是遥感最常用的波段; 3)红外线(IR):0.76-1000 μm。近红外0.76-3.0 μm’中红外3.0-6.0 μm; 远红外6.0-15.0 μm;超远红外15-1000 μm;(近红外又称光红外或反射红外;中红外和远红外又称热红外。) 4)微波:1mm-1m。全天候遥感;有主动与被动之分;具有穿透能力;发展潜力大。6.地物的反射光谱特性

ENVI遥感图像配准实验报告

ENVI遥感图像配准 一、实验目的: 1、掌握ENVI软件的基本操作和对图像进行基本处理,包括打开图像,保存图像。 2、初步了解图像配准的基本流程及采用不同校准及采样方法生成匹配影像的特点。 3、深刻理解和巩固基本理论知识,掌握基本技能和动手操作能力,提高综合分析问题的能力。 二、实验原理 (1)最邻近法 最邻近法是将最邻近的像元值赋予新像元。该方法优点是输出图像仍然保持原来图像的像元值,简单,处理速度快。缺点就是会产生半个像元位置偏移,可能造成输出图像中某些地物的不连贯。适用于表示分类或某种专题的离散数据,如土地利用,植被类型等。

双线性插方法是使用临近4个点的像元值,按照其距插点的距离赋予不同的权重,进行线性插。该方法具有平均化的滤波效果,边缘受到平滑作用,而产生一个比较连贯的输出图像,其缺点是破坏了原来的像元值,在后来的波谱识别分类分析中,会引起一些问题。 示意图: 由梯形计算公式: 故 同理 最终得:

三次卷积插法是一种精度较高的方法,通过增加参与计算的邻近像元的数目达到最佳的重采样结果。使用采样点到周围16邻域像元距离加权计算栅格值,方法与双线性插相似,先在Y 方向插四次(或X 方向),再在X 方向(或Y 方向)插四次,最终得到该像元的栅格值。该方法会加强栅格的细节表现,但是算法复杂,计算量大,同样会改变原来的栅格值,且有可能会超出输入栅格的值域围。适用于航片和遥感影像的重采样。 作为对双线性插法的改进,即“不仅考虑到四个直接邻点灰度值的影响,还考虑到各邻点间灰度值变化率的影响”,立方卷积法利用了待采样点周围更大邻域像素的灰度值作三次插值。其三次多项式表示为: 我们可以设需要计算点的灰度值f(x,y)为:

数字图像处理matlab代码

一、编写程序完成不同滤波器的图像频域降噪和边缘增强的算法并进行比较,得出结论。 1、不同滤波器的频域降噪 1.1 理想低通滤波器(ILPF) I1=imread('eight.tif'); %读取图像 I2=im2double(I1); I3=imnoise(I2,'gaussian',0.01); I4=imnoise(I3,'salt & pepper',0.01); figure,subplot(1,3,1); imshow(I2) %显示灰度图像 title('原始图像'); %为图像添加标题 subplot(1,3,2); imshow(I4) %加入混合躁声后显示图像 title('加噪后的图像'); s=fftshift(fft2(I4)); %将灰度图像的二维不连续Fourier 变换的零频率成分 移到频谱的中心 [M,N]=size(s); %分别返回s的行数到M中,列数到N中n1=floor(M/2); %对M/2进行取整 n2=floor(N/2); %对N/2进行取整 d0=40; %初始化d0 for i=1:M for j=1:N d=sqrt((i-n1)^2+(j-n2)^2); %点(i,j)到傅立叶变换中心的距离 if d<=d0 %点(i,j)在通带内的情况 h=1; %通带变换函数 else %点(i,j)在阻带内的情况 h=0; %阻带变换函数 end s(i,j)=h*s(i,j); %ILPF滤波后的频域表示

end end s=ifftshift(s); %对s进行反FFT移动 s=im2uint8(real(ifft2(s))); %对s进行二维反离散的Fourier变换后,取复 数的实部转化为无符号8位整数 subplot(1,3,3); %创建图形图像对象 imshow(s); %显示ILPF滤波后的图像 title('ILPF滤波后的图像(d=40)'); 运行结果: 1.2 二阶巴特沃斯低通滤波器(BLPF) I1=imread('eight.tif'); %读取图像 I2=im2double(I1); I3=imnoise(I2,'gaussian',0.01); I4=imnoise(I3,'salt & pepper',0.01); figure,subplot(1,3,1); imshow(I2) %显示灰度图像 title('原始图像'); %为图像添加标题 subplot(1,3,2); imshow(I4) %加入混合躁声后显示图像 title('加噪后的图像'); s=fftshift(fft2(I4));%将灰度图像的二维不连续Fourier 变换的零频率成分 移到频谱的中心 [M,N]=size(s); %分别返回s的行数到M中,列数到N中n=2; %对n赋初值

遥感数字图像处理教程实习报告

遥感数字图像处理教程实习报告

《数字图像处理》 课程实习报告 ( 2011 - 2012学年第 1 学期) 专业班级:地信09-1班 姓名:梁二鹏 学号:310905030114 指导老师:刘春国 ---------------------------------------------- 实习成绩: 教师评语: 教 师

签 名 : 年月日 实习一:图像彩色合成实习 一、实验目的 在学习遥感数字图像彩色合成基础上,应用所学知识,基于遥感图像处 理软件ENVI进行遥感数字图像彩色合成。 二、实验内容 彩色合成:利用TM图像can_tmr.img,实现灰度图像的密度分割、多波 段图像的真彩色合成、假彩色合成和标准假彩色合成。 三、实验步骤 1、显示灰度图像主要步骤: 1、打开ENVI4.7,单击FILE菜单,在下拉菜单中选择open image file 选 项,然后在弹出的对话框中选择can_tmr.img文件,单击打开。 2、在可用波段列表对话框中,选中某一波段图像,选中gray scale单选按 钮,单击LOAD BAND按钮,显示一幅灰度图像。 3、在可用波段列表对话框中,选择其他某一波段图像,进行显示。

4、利用可用波段列表中的display按钮,同时有多个窗口显示多个波段图像。 5、链接显示。利用图像窗口tool菜单下的link子菜单link display实现多图 像的链接显示。如图所示:红色方框。 6、使用tool菜单下的Cursor Location/value和pixel Locator功能在确定像 素的值和位置。

遥感地学分析读书报告

成像光谱技术研究动态 王立平刘洪博 1 引言 地物的反射辐射光谱特征是遥感的主要物理基础,是开展地球表层物质的物性和空间结构分析,进而加以识别的主要依据。成像光谱技术具有高光谱分辨率、超多波段和图谱合一的特点,在大尺度范围内探测地表物质连续光谱特性的同时,还获取了地物的空间形态和状态信息。成像光谱仪的光谱分辨率越高,所反映地物光谱特征就越精细,甚至可获取与实验室或地面实测光谱类似的曲线,为地物或地物成份的遥感识别奠定了基础。 2 成像光谱技术的发展与现状 成像光谱遥感所用的仪器是成像光谱仪。从世界范围来看,美国的成像技术发展较早,也最具代表性。从20世纪80年代到现在,美国已经研制了三代成像光谱仪。 第一代成像光谱仪的代表是航空成像光谱仪AIS。它由美国国家航空和航天管理局NASA所属的喷气推进实验室JPL设计,已于1984-1986年装在NASA的C-130飞机上飞行。这是一台装有二维、近红外阵列探测器的实验仪器,128个通道,光谱覆盖范围从1.2~2.4μm,并在内华达Cuprite地区的应用中取得很好的效果。 第二代成像光谱仪的代表是机载可见光/近红外成像光谱仪AVIRIS,它有224个通道,使用光谱范围为0.41~2.45μm,每个通道的波段宽约为10nm。曾放在改装后的高空U2飞机上使用.为目前最常用的航空光谱仪之一。 基于NASA仪器的成功应用,也基于采矿工业及石油工业的需求,在AVIRIS之后,地球物理环境研究公司GER又研制了l台64通道的高光谱分辨率扫描仪GERIS。其中63个通道为高光谱分辨率扫描仪,第64通道是用来存储航空陀螺信息。该仪器由3个单独的线性阵列探测器的光栅分光计组成。它与其他仪器的区别是在不同的光谱范围区内,通道的光谱宽度是不同的。

遥感图像实验报告

遥感图像实验报告 一.实验目的 1、初步了解目前主流的遥感图象处理软件ERDAS的主要功能模块。 2、掌握Landsat ETM遥感影像数据,数据获取手段.掌握遥感分类的方法, 土地利用变化的分析,植被变化分析,以及利用遥感软件建模的方法。 3、加深对遥感理论知识理解,掌握遥感处理技术平台和方法。 二.实验内容 1、遥感图像的分类 2、土地利用变化分析,植被变化分析 3、遥感空间建模技术 三.实验部分 1.遥感图像的分类 (1)类别定义:根据分类目的、影像数据自身的特征和分类区收集的信息确定分类系统; (2)特征判别:对影像进行特征判断,评价图像质量,决定是否需要进行影像增强等预处理; (3)样本选择:为了建立分类函数,需要对每一类别选取一定数目的样本;(4)分类器选择:根据分类的复杂度、精度需求等确定哪一种分类器; (5)影像分类:利用选择的分类器对影像数据进行分类,有的时候还需要进行分类后处理;分类图如下:

图1.1 1992年土地利用图 图1.2 2001年土地利用图

(6)结果验证:对分类结果进行评价,确定分类的精度和可靠性。 图1.3 1992年精度图 图1.4 2002年精度图 2.土地利用变化 2.1 两年土地利用相重合区域 (1)在两年的遥感影像中选择相同的区域。 Subset(x:568121~684371,y:3427359~3288369),过程如下:

图2.1 截图过程图 图2.2.2 截图过程图

(2)土地利用专题地图如下: 图2.2.3 1992年专题地图 图2.2.4 2001年土地利用图

图像处理实例(含Matlab代码)

信号与系统实验报告——图像处理 学院:信息科学与工程学院 专业:2014级通信工程 组长:** 组员:** 2017.01.02

目录 目录 (2) 实验一图像一的细胞计数 (3) 一、实验内容及步骤 (3) 二、Matlab程序代码 (3) 三、数据及结果 (4) 实验二图像二的图形结构提取 (5) 一、实验内容及步骤 (5) 二、Matlab程序代码 (5) 三、数据及结果 (6) 实验三图像三的图形结构提取 (7) 一、实验内容及步骤 (7) 二、Matlab程序代码 (7) 三、数据及结果 (8) 实验四图像四的傅里叶变化及巴特沃斯低通滤波 (9) 一、实验内容及步骤 (9) 二、Matlab程序代码 (9) 三、数据及结果 (10) 实验五图像五的空间域滤波与频域滤波 (11) 一、实验内容及步骤 (11) 二、Matlab程序代码 (11) 三、数据及结果 (12)

实验一图像一的细胞计数 一、实验内容及步骤 将该图形进行一系列处理,计算得到途中清晰可见细胞的个数。 首先,由于原图为RGB三色图像处理起来较为麻烦,所以转为灰度图,再进行二值化化为黑白图像,得到二值化图像之后进行中值滤波得到细胞分布的初步图像,为了方便计数对图像取反,这时进行一次计数,发现得到的个数远远多于实际个数,这时在进行一次中值滤波,去掉一些不清晰的像素点,剩下的应该为较为清晰的细胞个数,再次计数得到大致结果。 二、Matlab程序代码 clear;close all; Image = imread('1.jpg'); figure,imshow(Image),title('原图'); Image=rgb2gray(Image); figure,imshow(Image),title('灰度图'); Theshold = graythresh(Image); Image_BW = im2bw(Image,Theshold); Reverse_Image_BW22=~Image_BW; figure,imshow(Image_BW),title('二值化图像'); Image_BW_medfilt= medfilt2(Image_BW,[3 3]); figure,imshow(Image_BW_medfilt),title('中值滤波后的二值化图像'); Reverse_Image_BW = ~Image_BW_medfilt; figure,imshow(Reverse_Image_BW),title('图象取反'); Image_BW_medfilt2= medfilt2(Reverse_Image_BW,[20 20]); figure,imshow(Image_BW_medfilt2),title('第二次中值滤波的二值化图像'); [Label, Number]=bwlabel(Image_BW_medfilt,8);Number [Label, Number]=bwlabel(Image_BW_medfilt2,8);Number

遥感数字图像处理实验报告

实验一 遥感图像统计特性 一、实验目的 掌握遥感图像常用的统计特性的意义和作用,能运用高级程序设计语言实现遥感图像统 计参数的计算。 二、实验内容 编程实现对遥感图像进行统计特性分析,均值、方差(均方差)、直方图、相关系数等。 三、实验原理 1.均值 像素值的算术平均值,反映图像中地物的平均反射强度。 11 00 (,) N M j i f i j f MN --=== ∑∑ 2.方差(或标准差) 像素值与平均值差异的平方和,反映了像素值的离散程度。也是衡量图像信息量大小的 重要参数。 11 2 00 2[(,)] N M j i f i j f MN σ--==-= ∑∑ 3. 相关系数 反映了两个波段图像所包含信息的重叠程度。f , g 分别为两个波段的图像,它们之间的 相关系数计算公式为: 11 [((,))((,))] (,)M N f g f i j e g i j e C f g ---?-= ∑∑ 其中, e f , e g 分别为两个波段图像的均值。 四、实验步骤和内容 1.实验代码 clc clear all I =imread ('m1.jpg'); whos I %显示图像信息 figure (1),imshow (I ); R =double (I (:,:,1)); G =double (I (:,:,2)); B =double (I (:,:,3)); %求图像的R,G,B 的均值,avg=mean(mean(I))

%求图像的R,G,B的均值 mean(R(:)) mean(G(:)) mean(B(:)) %求R,G,B的方差 varR=var(R(:)); varG=var(G(:)) varB=var(B(:)) %求RG,RB,GB的相关系数 corrcoef(R(:),G(:)) corrcoef(R(:),B(:)) corrcoef(B(:),G(:)) 2.原始图像 Figure 1原始图像3.实验结果 R,G,B的均值

遥感图像预处理实验报告

实验前准备:遥感图像处理软件认识 1、实验目的与任务: ①熟悉ENVI软件,主要是对主菜单包含内容的熟悉; ②练习影像的打开、显示、保存;数据的显示,矢量的叠加等。 2、实验设备与数据 设备:遥感图像处理系统ENVI4.4软件; 数据:软件自带数据和河南焦作市影响数据。 3、实验内容与步骤: ⑴ENVA软件的认识 如上图所示,该软件共有12个菜单,每个菜单都附有下拉功能,里面分别包含了一些操作功能。 ⑵打开一幅遥感数据 选择File菜单下的第一个命令,通过该软件自带的数据打开遥感图像,可知,打开一幅遥感影像有两种显示方式。一种是灰度显示,另一种是RGB显示。 Gray(灰度显示)RGB显示 ⑶保存数据 ①选择图像显示上的File菜单进行保存; ②通过主菜单上的Save file as进行保存

⑷光谱库数据显示 选择Spectral > Spectral Libraries > Spectral Library Viewer。将出现Spectral Library Input File 对话框,允许选择一个波谱库进行浏览。点 击“Open Spectral Library”,选择某一所需的 波谱库。该波谱库将被导入到Spectral Library Input File 对话框中。点击一个波谱库的名称, 然后点击“OK”。将出现Spectral Library Viewer 对话框,供选择并绘制波谱库中的波谱曲线。 ⑸矢量化数据 点选显示菜单下的Tools工具栏,接着选择下面的第四个命令,之后选择第一个命令,对遥感图像进行矢量化。点击鼠标左键进行区域选择,选好之后双击鼠标右键,选中矢量化区域。 ⑹矢量数据与遥感影像的叠加与切割 选择显示菜单下的Tools工具,之后点选第一个 Link命令,再选择其下面的第一个命令,之后 OK,结束程序。 选择主菜单下的Basic Tools 菜单,之后选择 其中的第二个命令,在文件选择对话框中,选择 输入的文件(可以根据需要构建任意子集),将 出现Spatial Subset via ROI Parameters 对 话框通过点击矢量数据名,选择输入的矢量数 据。使用箭头切换按钮来选择是否遮蔽不包含在 矢量数据中的像元。 遥感图像的辐射定标 1、实验目的与任务: ①了解辐射定标的原理; ②使用ENVI软件自带的定标工具定标; ③学习使用波段运算进行辐射定标。 2、实验内容与步骤: ⑴辐射定标的原理 辐射定标就是将图像的数字量化值(DN)转化为辐射亮度值或者反射率或者表面温度等

基于MATLAB图像处理报告

基于M A T L A B图像处理报告一、设计题目 图片叠加。 二、设计要求 将一幅礼花图片和一幅夜景图片做叠加运算,使达到烟花夜景的美图效果。 三、设计方案 、设计思路 利用matlab强大的图像处理功能,通过编写程序,实现对两幅图片的像素进行线性运算,利用灰度变换的算法使图片达到预期的效果。 、软件介绍 MATLAB是matrix&laboratory两个词的组合,意为矩阵工厂(矩阵实验室)。是由美国mathworks公司发布的主要面对科学计算、可视化以及交互式程序设计的高科技计算环境。它将数值分析、矩阵计算、科学数据可视化以及非线性动态系统的建模和仿真等诸多强大功能集成在一个易于使用的视窗环境中,为科学研究、工程设计以及必须进行有效数值计算的众多科学领域提供了一种全面的解决方案,并在很大程度上摆脱了传统非交互式程序设计语言(如C、Fortran)的编辑模式,代表了当今国际科学计算软件的先进水平。 MATLAB和Mathematica、Maple并称为三大数学软件。它在数学类科技应用软件中在数值计算方面首屈一指。MATLAB可以进行矩阵运算、绘制函数和数据、实现算法、创建用户界面、连接其他编程语言的程序等,主要应用于工程计算、控制设计、信号处理与通讯、图像处理、信号检测、金融建模设计与分析等领域。 MATLAB的基本数据单位是矩阵,它的指令表达式与数学、工程中常用的形式十分相似,故用MATLAB来解算问题要比用C,FORTRAN等语言完成相同的事情简捷得多,并且MATLAB 也吸收了像Maple等软件的优点,使MATLAB成为一个强大的数学软件。在新的版本中也加入了对C,FORTRAN,C++,JAVA的支持。可以直接调用,用户也可以将自己编写的实用程序导入到MATLAB函数库中方便自己以后调用,此外许多的MATLAB爱好者都编写了一些经典的程序,用户直接进行下载就可以用。

遥感原理与应用实习

学号xxxx 天津城建大学 实习报告 遥感原理与应用实习 起止日期:2013 年12 月23日至2014年1月3 日 学生姓名Xx 班级XX 成绩 指导教师(签字) XX学院 2014年1 月3日

一、实习目的 “遥感原理与图像处理”课程是测绘工程专业的一门重要专业课,遥感信息是测绘、资源调查、环境监测、灾害评价诸方面应用的主要数据源。已在科学研究、工农业生产、军事、公安、医疗卫生、教育等许多领域得到广泛应用,产生了巨大的经济效益和社会效益,对推动社会发展,改善人们生活水平都起到了重要作用。未来要建立的数字地球是对真实地球及其相关现象数字化描述的一个虚拟地球。遥感技术将为数字地球提供动态的高分辨率、高光谱影像,用遥感影像生成的三维数字地面模型(DEM),以及地物和环境的各种属性数据等一些数字地球中最基础的数据。 “遥感实习”目的是培养学生进行遥感技术应用和图像数字处理的实际操作能力。要求了解一些基本的地物波谱反射率的野外测定方法,理解遥感图像目视解译,了解航天(或航空)像片识读与野外调绘。 二、实验项目基本要求 1.熟悉一种遥感图像处理软件 2.遥感影像的认知,进行图像剪切,波段组合与图像显示 3.图像的几何校正 4.遥感影像增强处理 5.遥感影像解译 三、实习步骤(包括原理,方法,操作过程) 1.图象剪切, 波段组合与图像显示 原图像比较大,数据量大处理不方便,对齐剪切便于计算机处理,也能达到实习目的 剪切DatePrep>SubsetImage命令如下图所示

波段组合Raster>Band Combinations 打开波段设置对话框 1)真彩色合成,即3、2、1波段分别赋予红、绿、蓝色,则获得自然彩色合成图 像,图像的色彩与原地区或景物的实际色彩一致.如下图 2)标准假彩色合成,即4、3、2波段分别赋予红、绿、蓝色,获得图像植被成红 色,由于突出表现了植被的特征,应用十分的广泛,而被称为标准假彩色。

遥感数字图像处理实习报告

目录 一、实习目的 (1) 二、实习要求 (1) 三、实习内容 (1) 3.1 图像裁剪 (1) 3.2 图像配准 (3) 3.3 监督分类 (6) 四、实习总结 (16)

一、实习目的 (1)熟悉使用erdas软件 (2)了解图像处理的原理与步骤,提高动手能力 二、实习要求 (1)对图像进行裁剪 (2)将TM图像与与临川区图像进行配准(校正) (3)对图像进行监督分类 三、实习内容 3.1 图像裁剪 利用临川区行政边界AOI文件对TM图象进行裁剪,裁剪出临川区TM图象。 基本步骤如下: (1)打开erdas imagine 9.2,单击“data proparation”,如图1 图1 data proparation 界面 (2)单击“subset imagine”,进入以下界面: 图2 subset界面

(3)图像裁剪的方法有三种 ①在TM影像上,通过AOI工具进行裁剪,裁剪出所需要的范围 ②通过使用文件裁剪 ③ AOI裁剪,在subset界面,单击AOI, 图3 AOI文件的输入 (4)输入AOI文件后,确认并进行裁剪,得到以下结果 图4 裁剪后的结果

3.2 图像配准 map-to-image: 1:10万临川区土地利用图与TM图象配准; 要求最初选GCP点6-10个,及检测点5个,各点均匀分布,RMS检验误差小于30米(1个像元)。 (1)进入“data proparation”界面,单击以下图中所示: 图5 data proparation 界面 (2)选择视窗,应该注意影像图的选择,选择tiff格式的临川区规划图(需校正的影像) 图6 视窗的选择 (3)模型的选择(多项式) 图7 模型选择界面

相关主题