档案馆借的地形数据是分幅的,使用python批量镶嵌
使用 arcpy.MosaicToNewRaster_management
方法。
python
# -*- coding: UTF-8 -*-
import arcpy
import glob
import os
'''
数据目录如下:
dem
|-- 分幅1
tif数据
|-- 分幅2
tif数据
'''
'''
思路:
两个输入路径 -> 两个范围的栅格数据
'''
arcpy.CheckOutExtension("Spatial") # 检查是否注册该模块
# 输入路径 应该注意,中文路径,会导致读不出文件
inws = r"D:\dem\G48"
# 输出路径
outws = r"Z:\lm\dem"
# 以下逻辑为:
# 获得子目录全路径,并进入子目录
path_list = os.listdir(inws)
rasters = []
# 进入二级目录
for i in range(len(path_list)):
path = inws + "\\" + path_list[i]
ras= glob.glob(os.path.join(path, "*.bil"));
rasters.append(ras[0])
# 因为不太熟悉python,在for循环里没法修改全局变量
# 4090是指该文件夹下的全部数据
if(i == 4090):
base = rasters[0] # 指定该工作空间下的一副影像为基础影像,为后面的参数提取做准备
out_coor_system = arcpy.Describe(base).spatialReference # 坐标系统
dataType = arcpy.Describe(base).DataType # 数据类型
piexl_type = arcpy.Describe(base).pixelType # 像素类型
cellwidth = arcpy.Describe(base).meanCellWidth # 像元宽度
bandcount = arcpy.Describe(base).bandCount # 获取bandCount
mosaic_ras = [] # 待拼接影像的文件名,且中间以;隔开,例如:a.tif;b.tif;c.tif
for ras in rasters:
mosaic_ras.append(ras)
mosaic_ras = ";".join(rasters) # 字符串拼接
nameT = 'hebing' + '.bil'
# 执行拼接操作
arcpy.MosaicToNewRaster_management(mosaic_ras, outws, nameT, out_coor_system,
"32_BIT_FLOAT",
cellwidth, bandcount, "LAST", "FIRST")
print os.path.basename(nameT) + " ---- OK! ---- "
print(" --- All project is OK! --- ")