本篇文章將說明如何利用ArcGIS 10.1自帶的Python IDLE進(jìn)行遙感影像的批量拼接與裁剪。

1.運(yùn)行環(huán)境：ArcGIS10.1 (安裝傳送門)、Python IDLE

2.數(shù)據(jù)來源：地理空間數(shù)據(jù)云 GDEMV2 30M分辨率數(shù)字高程數(shù)據(jù)

3.解決問題：制作山西省的DEM影像

如下圖所示，以30M分辨率數(shù)字高程數(shù)據(jù)為例，影像皆是固定范圍的經(jīng)緯度保存在其服務(wù)器上，外在表現(xiàn)以小幅正方形影像。如果手動(dòng)進(jìn)行拼接，工作量會(huì)非常大且容易出錯(cuò)。

我們的目標(biāo)：批處理，寫一次代碼，處理多幅影像

1.查找目標(biāo)范圍的經(jīng)緯度信息，本文以山西為例，經(jīng)緯度范圍在N33-N41, E109-E114直接，所以先將這個(gè)范圍內(nèi)的影像下載后，解壓后放在一個(gè)文件夾下。

2.拼接影像

開始菜單-ArcGIS-IDLE（Python GUI） 打開IDLE。File-New Window, 文件以 .py結(jié)尾，如MosaicToNewRasters.py

?

import arcpy import os#指定工作目錄，即存放影像的目錄 arcpy.env.workspace = r"E:\arcpyData\shanxi\N34_N35"#指定該工作空間下的一副影像為基礎(chǔ)影像，為后面的參數(shù)提取做準(zhǔn)備 base = "ASTGTM2_N34E109_dem.tif"#以下一段代碼是為執(zhí)行拼接做參數(shù)準(zhǔn)備 out_coor_system = arcpy.Describe(base).spatialReference #獲取坐標(biāo)系統(tǒng) dataType = arcpy.Describe(base).DataType piexl_type = arcpy.Describe(base).pixelType cellwidth = arcpy.Describe(base).meanCellWidth #獲取柵格單元的的寬度 bandcount = arcpy.Describe(base).bandCount #獲取bandCount#打印一些信息 print out_coor_system.name print dataType print piexl_type print cellwidth print bandcountarcpy.CheckOutExtension("Spatial")#提取待拼接影像的文件名，且中間以;隔開，例如：a.tif;b.tif;c.tif rasters = [] for ras in arcpy.ListRasters("*dem.tif"): #for循環(huán)，對(duì)wrokspace下的所有以dem.tif結(jié)尾的影像進(jìn)行過濾rasters.append(ras)ras_list = ";".join(rasters) #字符串拼接#打印出來，看看什么結(jié)果吧 print ras_list#指定輸出文件夾 outFolder = r"E:\arcpyData\shanxi\N34_N35\n34_n35_out"#執(zhí)行拼接操作 arcpy.MosaicToNewRaster_management(ras_list, outFolder, "shanxi_n34_n35_dem_data.tif", out_coor_system, "16_BIT_SIGNED", cellwidth, bandcount, "LAST", "FIRST")

RUN-->Run Module F5 執(zhí)行代碼

?

這里用到了arcpy包下面的 MosaicToNewRaster_management()函數(shù)

得到如下結(jié)果

3.以山西省shp邊界裁剪影像

再新建一個(gè)Window，文件命為BatchExtractByMask.py

?

import arcpy import glob import osarcpy.CheckOutExtension('Spatial')#指定先前拼接后的遙感影像所在目錄 inws = r"E:\arcpyData\shanxi\shanxi_regular_dem\shanxi_province_dem"#指定裁剪后的影響存放目錄 outws = r"E:\arcpyData\shanxi\shanxi_regular_dem\shanxi"#指定shp范圍邊界文件，即目標(biāo)區(qū)域的邊界 mask = r"E:\arcpyData\shanxi\shanxi_shp\shanxi.shp"#利用glob包，將inws下的所有tif文件讀存放到rasters中 rasters = glob.glob(os.path.join(inws, "*.tif"))#循環(huán)rasters中的所有影像，進(jìn)行按掩模提取操作 for ras in rasters:outname = os.path.join(outws, os.path.basename(ras).split(".")[0] + "_clp.tif") #指定輸出文件的命名方式（以被裁剪文件名+_clip.tif命名）out_extract = arcpy.sa.ExtractByMask(ras, mask) #執(zhí)行按掩模提取操作out_extract.save(outname) #保存數(shù)據(jù)

?

執(zhí)行腳本代碼，裁剪后的影像被存放在outw所指定的文件夾

4.通過ArcMap加載影像，最終結(jié)果：

這就實(shí)現(xiàn)了批量拼接與裁剪影像的工作，其實(shí)拼接與裁剪無所謂先后順序。

?

如果你覺得本文對(duì)你有幫助，是支持贊賞的哦：）

?

如遇到問題，歡迎通過公眾號(hào)留言給作者，以便共同探討。

郵箱：thinkingingis@qq.com

微信公眾號(hào)：

總結(jié)

以上是生活随笔為你收集整理的利用ArcGIS Python批量拼接裁剪遥感影像(arcpy batch processing)的全部?jī)?nèi)容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網(wǎng)站內(nèi)容還不錯(cuò)，歡迎將生活随笔推薦給好友。