??本文介紹基于Python語(yǔ)言arcpy模塊,實(shí)現(xiàn)柵格影像圖層建立與多幅遙感影像數(shù)據(jù)批量拼接(Mosaic)的操作。
??首先,相關(guān)操作所需具體代碼如下:
import os
import arcpy
file_path="G:/Postgraduate/LAI_Glass_RTlab/A2018161_Dif/DRT/"
out_file_path="G:/Postgraduate/LAI_Glass_RTlab/A2018161_Dif/DRT/"
out_file_name="Global.tif"
file_name_list=os.listdir(file_path)
tif_file_path=file_path+file_name_list[0]
cell_size_x=arcpy.GetRasterProperties_management(tif_file_path,"CELLSIZEX")
cell_size=cell_size_x.getOutput(0)
value_type=arcpy.GetRasterProperties_management(tif_file_path,"VALUETYPE")
describe=arcpy.Describe(tif_file_path)
spatial_reference=describe.spatialReference
arcpy.CreateRasterDataset_management(out_file_path,out_file_name,cell_size,"16_BIT_SIGNED",
spatial_reference,"1")
out_file=out_file_path+out_file_name
for file in file_name_list:
file_path_name=file_path+file
print(file_path_name)
arcpy.Mosaic_management([file_path_name],out_file)
??其中,file_path為存放有多景初始遙感影像的路徑格式為.tif柵格文件(如果不是.tif格式,例如是.hdf等文件,需首先進(jìn)行文件格式的轉(zhuǎn)換);out_file_path為拼接后所得結(jié)果柵格圖層的存放路徑;out_file_name為拼接后所得結(jié)果柵格圖層的文件名稱,其可選格式有很多,如下圖所示。

??在這里,我們默認(rèn)所得拼接結(jié)果圖層為一個(gè)(也就是file_path文件夾中全部的待處理遙感影像最終全拼接在一起);如果大家需要使得拼接結(jié)果圖層是多幅(也就是file_path文件夾中待處理遙感影像依據(jù)區(qū)域、時(shí)間等分為很多不同的部分,每一部分拼接在一起),可以參考這篇博客(https://blog.csdn.net/zhebushibiaoshifu/article/details/118878435),利用其中的循環(huán)方式實(shí)現(xiàn)需求。
??隨后,通過os.listdir()函數(shù)獲取file_path路徑下的柵格文件,并存儲(chǔ)于file_name_list列表中。
??接下來(lái)需要?jiǎng)?chuàng)建一個(gè)新的柵格圖層。之所以要進(jìn)行這一步驟,是因?yàn)楸疚暮笃谶x擇用arcpy.Mosaic_management()函數(shù)進(jìn)行柵格的批量拼接,因此需要首先創(chuàng)建一個(gè)新的、空的柵格圖層作為拼接的基準(zhǔn)。如果大家的需求不是批量拼接?xùn)鸥駭?shù)據(jù),而是單純想利用arcpy進(jìn)行新柵格的創(chuàng)建,那就只看這一部分的代碼即可。
??在這里,我們選擇用file_path路徑下的第一個(gè)柵格數(shù)據(jù)(下稱“第一柵格”)作為新柵格圖層中各項(xiàng)屬性(例如像素邊長(zhǎng)、像素?cái)?shù)據(jù)格式等)的依據(jù)。首先,arcpy.GetRasterProperties_management()函數(shù)獲取第一柵格的像素x邊邊長(zhǎng);因?yàn)橐话銝鸥駭?shù)據(jù)中像素都是正方形,因此我們就通過cell_size=cell_size_x.getOutput(0)將第一柵格的像素x邊邊長(zhǎng)作為新柵格圖層像素x邊與y邊二者的邊長(zhǎng)。再利用arcpy.GetRasterProperties_management()函數(shù)獲取第一柵格的數(shù)據(jù)格式;最后利用中間變量describe獲取第一柵格的空間參考信息。
??完成以上步驟后,將已獲取的第一柵格的各類信息通過函數(shù)arcpy.CreateRasterDataset_management()帶入新柵格中。在這里需要注意:盡可能在將要拼接時(shí)選擇新柵格為"16_BIT_SIGNED"及以下的數(shù)據(jù)格式(具體數(shù)據(jù)格式類別如下圖),且將file_path路徑下待拼接的柵格數(shù)據(jù)的數(shù)據(jù)格式也全部修改為這一格式;否則可能會(huì)由于數(shù)據(jù)量大而導(dǎo)致拼接過程極慢。我之前就是由于選用了32 bit float格式的柵格數(shù)據(jù)進(jìn)行拼接,導(dǎo)致全球范圍的MODIS一個(gè)植被產(chǎn)品數(shù)據(jù)拼接花了將近一天的時(shí)間。如果大家的柵格像素?cái)?shù)據(jù)包含小數(shù),可以通過乘上一個(gè)縮放系數(shù)的方式進(jìn)行數(shù)據(jù)整數(shù)化。

??代碼最后的一個(gè)for循環(huán),就是遍歷file_name_list中的各個(gè)柵格數(shù)據(jù),并通過arcpy.Mosaic_management()函數(shù)加以拼接即可。
??以上,便完成了本次批量拼接的操作。這里還有一點(diǎn)需要注意:由于arcpy模塊的限制,如果大家的Python版本是3.0及以上,往往不能直接運(yùn)行上述代碼,最好是在ArcMap的Python運(yùn)行框或其對(duì)應(yīng)IDLE(如下圖所示)中運(yùn)行。
