遥感影像的镶嵌与融合 PPT
- 格式:ppt
- 大小:394.00 KB
- 文档页数:15
遥感实验五数字图像处理------------图像镶嵌、裁切及融合一、实验目的学会图像镶嵌、图像裁切及图像融合等技术,通过实际影像的操作,制作可用于实际工作的某区域遥感图像,为下一次实验准备数据。
二、实验数据某区域的遥感图像:11942E20000504.rar、11943E20010304.rar;某区域的范围:xianyou.shp三、实验内容及主要步骤1、图像镶嵌注:要镶嵌的两幅或多幅影像要求具有相同的投影信息,如果不同,则需要首先统一。
ERDAS IMAGINE中提供了投影转换的工具,点击、选择Reproject Images;或者,也可以在ArcGIS的ArcToolbox中选择Projections and Transformations/Raster/Project Raster进行转换。
以下以ERDAS IMAGINE 软件为例进行投影转换。
1.1.投影定义和转换在ERDAS中,点击DataPrep,在下拉选项卡中点击Rejection Images,在Input File中输入需要进行投影转换的影像数据——福建某地区2000年5月30米分辨率的的多光谱影像(本例以TM4、3、2波段为例)。
在Output File设置保存路径和输出文件名。
在Categories中点击右侧的小地球标志进行投影定义。
投影参数设置如下图1.1示,点击OK,完成投影转换。
本图及以下各图均将WGS-84投影转换成Gauss Kruger投影。
图1.1同理,对裁切的多光谱小图进行投影转换,原理及步骤亦同上,图1.2示。
图1.2对全色波段影像数据tm11942_8进行投影变换,原理同多光谱影像投影变换,但在erdas 中进行投影转化时由于在选择categories时,选择了南半球国家投影类别发生错误,结果显示为一“倒像”,故tm11942_8影像采用ArcGIS软件进行投影转换,转换目的主要是讲投影信息中的Datum转成Krasovskv。
实验四遥感图像的拼接、裁剪、融合一、实习目的与要求·掌握图像拼接的原理,以及两幅图像拼接的时候需要的条件,掌握拼接技术;·学习通过ERDAS进行遥感图像规则分幅裁剪,不规则分幅裁剪的实验过程,能够对一幅大的遥感图像按照要求裁剪图像;·掌握不同分辨率图像的特性,详细理解各种融合方法的原理,以及各种融合方法的优缺点,能够根据不同的应用目的合理选择融合方法,掌握融合的操作过程;二、实验原理·图像拼接(mosaic image)是具有地理参考的若干相邻的图像合并成一幅图像或一组图像,需要拼接的图像必须含有地图投影也就是说图像必须经过几何校正处理,虽然所有的输入图像可以具有不同的投影类型,不同的象元大小,但必须有相同的波段数。
在进行图像拼接时需要确定一幅参考影像,参考图像作为图像拼接的基准,决定输出图像的地图投影和象元大小和数据类型。
·在实际工作中,经常需要根据研究区域的工作范围对图像进行分幅裁剪,erdas中可以对图像进行规则分幅裁剪(rectangle subset)和不规则分幅裁剪(pdygon subset),根据实际的应用对图像选择不同的裁剪方式。
·分辨率融合是对不同分辨率的摇杆图像进行融合处理,使处理后的图像既具有较好的空间分辨率又具有多光谱特征,从而增加图像的可解译性。
图像分辨率融合的关键是融合前两幅图像的配准以及融合方法的选择只有将不同空间分辨率的图像进行精确的配准才能达到满意的融合效果,而融合的方法的选择主要是由被融合图像的特性以及融合的目的进行选择的,同时需要对融合的原理有正确的认识。
三、实验内容和实验过程本次试验主要包括遥感图像拼接、遥感图像分幅裁剪、遥感图像分辨率融合。
下面分别介绍:1.图像拼接实验步骤:(1)启动图象拼接工具,在ERDAS图标面板工具条中,点击Dataprep/Data preparation/Mosaicc lmages—打开Mosaic Tool 视窗。
遥感图像的拼接和镶嵌1 ⾃定义镶嵌函数遥感图像的镶嵌,主要分为5⼤步骤:step1: 1)对于每⼀幅图像,计算其⾏与列;2)获取左上⾓X,Y3)获取像素宽和像素⾼4)计算max X 和 min Y,切记像素⾼是负值maxX1 = minX1 + (cols1 * pixelWidth)minY1 = maxY1 + (rows1 * pixelHeight)step2 :计算输出图像的min X ,max X,min Y,max YminX = min(minX1, minX2, …)maxX = max(maxX1, maxX2, …)y坐标同理step3:计算输出图像的⾏与列cols = int((maxX – minX) / pixelWidth)rows = int((maxY – minY) / abs(pixelHeight)step 4:创建⼀个输出图像driver.create()step 5:1)计算每幅图像左上⾓坐标在新图像的偏移值2)依次读⼊每幅图像的数据并利⽤1)计算的偏移值将其写⼊新图像中import osimport globimport mathimport gdaldef get_extent(fn):ds = gdal.Open(fn)gt = ds.GetGeoTransform()r = ds.RasterYSizec = ds.RasterXSizereturn (gt[0], gt[3], gt[0] + gt[1] * c, gt[3] + gt[5] * r)def mosiac(in_files):min_X, max_Y, max_X, min_Y = get_extent(in_files[0])for fn in in_files[1:]:minx, maxy, maxx, miny = get_extent(fn)min_X = min(min_X, minx)max_Y = max(max_Y, maxy)max_X = max(max_X, maxx)min_Y = min(min_Y, miny)ds1 = gdal.Open(in_files[0])band1 = ds1.GetRasterBand(1)transform1 = ds1.GetGeoTransform()pixelWidth1 = transform1[1]pixelHeight1 = transform1[5]bands = ds1.RasterCount# 获取输出图像的⾏与列cols = int((max_X - min_X) / pixelWidth1)rows = int((max_Y - min_Y) / abs(pixelHeight1))driver = ds1.GetDriver()dsOut = driver.Create(r'F:\algorithm\算法练习\拼接与镶嵌\mosiac1.tif', cols, rows, bands, band1.DataType)#1是bands,默认for file in in_files:ds2 = gdal.Open(file)rows2 = ds2.RasterYSizecols2 = ds2.RasterXSizedata2 = ds2.ReadAsArray(0, 0, cols2, rows2)transform2 = ds2.GetGeoTransform()minX2 = transform2[0]maxY2 = transform2[3]# 计算每张图左上⾓的偏移值(在输出图像中)xOffset2 = int((minX2 - min_X) / pixelWidth1)yOffset2 = int((maxY2 - max_Y) / pixelHeight1)for i in range( bands):dsOut.GetRasterBand(i+1).WriteArray(data2[i],xOffset2, yOffset2)geotransform = [min_X, pixelWidth1, 0, max_Y, 0, pixelHeight1]dsOut.SetGeoTransform(geotransform)dsOut.SetProjection(ds1.GetProjection())if__name__ == '__main__':os.chdir(r'F:\algorithm\算法练习\拼接与镶嵌\test2_next')in_files = glob.glob('*.tif')print(in_files)mosiac(in_files )2 采⽤gdal.Warp()提供的接⼝进⾏镶嵌def RasterMosaic():print("图像拼接")inputrasfile1 = gdal.Open(inputfilePath, gdal.GA_ReadOnly) # 第⼀幅影像inputProj1 = inputrasfile1.GetProjection()inputrasfile2 = gdal.Open(referencefilefilePath, gdal.GA_ReadOnly) # 第⼆幅影像inputProj2 = inputrasfile2.GetProjection()outputfilePath = 'G:/studyprojects/gdal/GdalStudy/Files/images/RasterMosaic2.tif'options=gdal.WarpOptions(srcSRS=inputProj1, dstSRS=inputProj1,format='GTiff',resampleAlg=gdalconst.GRA_Bilinear) gdal.Warp(outputfilePath,[inputrasfile1,inputrasfile2],options=options)。