需求说明
在进行DEM数据以及DSM数据生产时，可能会出现DEM数据内部部分点高程数据大于DSM数据的情况，对此，通常采用的做法是将DEM内这些点的高程数据用对应的DSM数据进行覆盖，在ArcGIS中，使用栅格计算器完成这一操作，栅格计算器的公式模板如下(将DEM和DSM替换为自己项目中的对应的数据)：
Con(DEM > DSM, DEM, DSM)
使用上述方法在面对多个数据（数十个DEM和DSM数据）时不太方便，为此，开发相应脚本以批量处理该操作
代码实现
核心函数
def dem_negative_calculator(in_dem,in_dsm,out_folder):
"""
去除DEM负值数据
:param in_dem:str，输入DEM数据文件夹
:param in_dsm:str，输入DSM数据文件夹
:param out_folder:str，输出结果文件夹
:return:
"""
arcpy.CheckOutExtension('Spatial')
dem_data_paths=[]
dsm_data_paths=[]
# 获取文件夹内所有DEM数据，项目中通常为 图幅号+"DEMK"的形式
for root, dirs, files in os.walk(in_dem):
for file in files:
if file.endswith('.img') and "DEM" in file.upper():
dem_data_paths.append(os.path.join(root, file))
dem_data_paths.sort()
# 获取文件夹内所有DSM数据，项目中通常为 图幅号+"DSMK"的形式
for root, dirs, files in os.walk(in_dsm):
for file in files:
if file.endswith('.img') and "DSM" in file.upper():
dsm_data_paths.append(os.path.join(root, file))
dsm_data_paths.sort()
# 设置进度条
arcpy.SetProgressor("step","Remove DEM negative value",0,len(dem_data_paths),1)
for dem_file in dem_data_paths:
dem_name=os.path.basename(dem_file)
# 更新进度条
arcpy.SetProgressorLabel("Processing {0}".format(dem_name))
arcpy.SetProgressorPosition()
# 在DSM数据中找到相应的DEM图幅号的数据
for dsm_file in dsm_data_paths:
dsm_name=os.path.basename(dsm_file)
# DEM和DSM图幅号数据的对比
if os.path.splitext(dem_name)[0][0:-4]==os.path.splitext(dsm_name)[0][0:-4]:
result_path=os.path.join(out_folder,dem_name)
# 使用栅格计算条件分析功能Con，公式模板：Con(DEM > DSM, DEM, DSM)
raster=arcpy.sa.Con(dem_file>dsm_file,dem_file,dsm_file)
if arcpy.Exists(result_path):
arcpy.Delete_management(result_path)
# 结果保存新的路径中
raster.save(result_path)
break
arcpy.ResetProgressor()
参数获取
dem_folder=arcpy.GetParameterAsText(0)
dsm_folder=arcpy.GetParameterAsText(1)
res_folder=arcpy.GetParameterAsText(2)
ArcGIS工具导入
参数设置：
运行界面：
