基于Acrgis的Python脚本在森林资源年度更新中的应用

时启龙，邱 琳，喻 俊

(江西省林业调查规划研究院，南昌 330046)

森林资源年度更新是摸清森林资源家底、进行森林资源管理的一项基础性工作[1]，也是掌握林地及林木资源空间分布及动态变化情况，保持林地及林木资源数据准确性和时效性的有效手段，可为本地区生态文明试验区建设目标评价考核、科学发展观综合考核、领导干部自然资源资产离任审计、自然资源资产负债表编制等工作提供及时准确的基础数据，也为当地建设项目使用林地行政许可和林地保护行政执法、生态公益林与天然林资源保护管理等提供执法依据，并可为当地政府和林业主管部门科学决策、规范管理提供重要支撑。

森林资源年度更新每年都需要进行大量的地理空间数据处理工作，而且这种处理过程往往具有批量化和重复性的特点，在以往的数据处理中，一般是运用Arcgis软件对单一数据集或单一过程进行处理，不仅效率低下，而且容易出错[2]。Python是一种不受局限、跨平台的开源编程语言，它功能强大且简洁易读，可伸缩性强，既适用于大型项目，又可应用于一次性程序[3]，目前已延伸到Arcgis软件中，成为一种进行数据提取、数据转换、数据处理以及数据分析和地图制图自动化的脚本语言[4]，在Arcgis中运用Python语言可实现对地理数据处理的自动化和批量化，从而大大提高了工作效率，同时计算机视觉库OpenCV、三维可视化库VTK和医学图像处理库ITK等众多开源的科学计算软件包都提供了Python调用接口[5]。因此，我们在森林资源年度更新多个环节中进行了基于Acrgis的Python脚本批处理操作，也取得了满意效果。

1 遥感影像裁剪

森林资源年度更新一个重要的基础性工作就是以县域为单元对遥感影像进行裁剪，在以往的工作过程中，影像裁剪工作均是手动操作，首先需要参照影像图幅框中各影像的空间位置，在全省的影像中找到与某县域相交的所有影像，然后在Arcgis中进行裁剪，一个单元裁剪完成后，再手动进行下一个单元的裁剪配置，不仅工作量大，效率低，而且在查找影像过程中容易混淆出错。而采用Python脚本语言，则只需要设置好全省县界面、影像图幅框和全省遥感影像图的存放位置，就可以实现对全省一百多个单位和几百景遥感影像进行自动裁剪。具体实现代码如下：

#导入站点包和环境参数

import arcpy

from arcpy import env

import os

#定义初始条件

arcpy.CheckOutExtension("spatial")

arcpy.gp.overwriteOutput = 1

env.workspace = "c:/zygx"

xjm = "c:/zygx/xjm.shp"

tfk = "c:/zygx/tufukuang.shp"

xnum = "XIAN"

#加载全省县界面

arcpy.MakeFeatureLayer_management(xjm, "xjm")

#定义查询游标，遍历县界面中的元素

xcursors = arcpy.da.SearchCursor("xjm", ["XIAN","XIAN_NAME"])

for xcursor in xcursors:

qry= xnum + "=" + "'" + xcursor[0] + "'"

#选中裁剪单元

arcpy.SelectLayerByAttribute_management("xjm","NEW_SELECTION",qry)

#定义裁剪后影像的存放路径

folder = "c:/zygx/cj/ " + xcursor[1]

os.makedirs(folder)

#定义查询游标，遍历图幅框中的元素

qryy = "XIAN" + "= '" + xcursor[0] + "'"

tcursors = arcpy.da.SearchCursor(tfk, ["XIAN","TNAME"],qryy)

for tcursor in tcursors :

inpt = "c:/zygx/yx/" + tcursor[1] + ".ers"

outpt = folder + "/" + tcursor[1] + "_clip.img"

#裁剪影像

arcpy.Clip_management(inpt,"",outpt,"xjm","","NONE")

print(tcursor[1] + "Has Done")

print(xcursor[1] + "All Done")

print("All Done")

2 遥感影像拼接

遥感影像拼接作为影像处理过程中的一个重要环节，比较费事费力，而借助Arcgis的Python脚本语言可以实现对全省一百多个单位的遥感影像图进行自动化拼接，拼接过程简单易行，只需要设置好输入和输出数据的存放路径即可。具体实现代码如下：

import arcpy

from arcpy import env

import os

arcpy.CheckOutExtension("spatial")

arcpy.gp.overwriteOutput = 1

env.workspace = "c:/zygx/"

#遍历影像存放文件夹列表

filelists = os.listdir("c:/zygx/cj")

for filelist in filelists:

#获取文件夹中影像列表

env.workspace = "c:/zygx/cj/" + filelist

rasters = arcpy.ListRasters("*","img")

#影像拼接

yxNmae = filelist + ".img"

arcpy.MosaicToNewRaster_management(rasters, env.workspace, yxName, "", "8_BIT_UNSIGNED", "2", "3", "LAST", "FIRST")

3 森林资源小班数据库分发

森林资源小班数据库分发是森林资源年度更新的一项基础性准备工作，整个森林资源年度更新过程均在分发的森林资源小班数据库上进行。在以往的小班数据分发过程中，一般是先建立以单位代码和名称为文件名的文件地理数据库或个人地理数据库，然后在全省数据库中按照单位代码或单位名称选中本单位数据并导入所建地理数据库中，过程虽简单，但比较费时，且效率低。而借助Python脚本语言可以实现对这一操作的自动化运行。具体代码如下：

import arcpy

from arcpy import env

arcpy.CheckOutExtension("spatial")

arcpy.gp.overwriteOutput = 1

env.workspace = "c:/zygx/"

zym = "c:/zygx/jxzym2018.gdb/jxzym2018"

xjm = "c:/zygx/xjm.shp"

zyname = u'资源面2018'

xnum = "XIAN"

#加载森林资源小班数据库

arcpy.MakeFeatureLayer_management(zym,"zym")

#将加载的资源小班数据库隐藏不显示

mxd = arcpy.mapping.MapDocument("CURRENT")

ryx = arcpy.mapping.ListDataFrames(mxd,"layers")[0]

for lyr in arcpy.mapping.ListLayers(mxd, "", ryx):

if lyr.name == "zym":

lyr.visible = "false"

#定义查询游标，遍历数据库中元素

cursors = arcpy.da.SearchCursor(xjm,["XIAN","XIAN_NAME"])

for cursor in cursors:

#建立文件地理数据库

arcpy.CreateFileGDB_management("c:/zygx/zysj/",cursor[0] + cursor[1])

#选中森林资源小班数据库

qry = xnum + "=" + "'" + cursor[0] + "'"

arcpy.SelectLayerByAttribute_management("zym","NEW_SELECTION",qry)

#导出所选数据

outpt = "c:/zygx/zysj/" + cursor[0] + cursor[1] + ".gdb"

tname = cursor[1] + zyname

arcpy.FeatureClassToFeatureClass_conversion("zym",outpt,tname)

print(cursor[1] + "Has Done")

print("All Done")

4 外业调查图批量制图

调查图制作是森林资源调查中的一项重要工作，几乎应用于林业管理的各个方面，在传统的森林资源调查过程中，外业调查图纸的制作均是在设置好模板的基础上，一张一张的导出，而无法实现批量出图，这对于大批量的制图工作是不利的，不仅工作量巨大，而且需要占用大量的时间，时间成本较高。Python脚本语言为Arcgis定制了专门的批量地图制图功能，不仅可以利用已有的模板制图，而且可以重新设定制图模板。由于森林资源年度更新外业调查图图式一致，本项目则借助已有的制图模板，利用Python脚本的DataDrivenPages类并结合Arcgis中的Data Driven Pages工具，来实现对外业调查图纸的批量化和自动化制作。具体实现代码如下：

import arcpy

from arcpy import env

arcpy.CheckOutExtension("spatial")

arcpy.gp.overwriteOutput = 1

env.workspace = "c:/zygx/"

#定义当前文档

mxd = arcpy.mapping.MapDocument("CURRENT")

#依据标准图幅框切换制图范围

for pageNum in range(1,mxd.dataDrivenPages.pageCount):

mxd.dataDrivenPages.currentPageID=pageNum

mapName = mxd.dataDrivenPages.pageRow.getValue(mxd.dataDrivenPages.pageNameField.name)

raster = "c:/zygx/ditu" + mapName + ".tif"

#依据图幅框加载底图

arcpy.MakeRasterLayer_management(raster,mapName)

#修改图名

elm = arcpy.mapping.ListLayoutptElements(mxd)[0]

for elm in arcpy.mapping.ListLayoutptElements(mxd):

if elm.name == "ZTH":

elm.text = mapName

#导出地图

arcpy.mapping.ExportToJPEG(mxd, r"c:/zygx/zhitu/" + mapName + ".jpg",resolution = 300)

#保存MXD文件

scy = r"c:/zygx/zhitu/" + mapName + ".mxd"

mxd.saveACopy(scy)

#删除所加载底图

ryx = arcpy.mapping.ListDataFrames(mxd,"layers")[0]

for lyr in arcpy.mapping.ListLayers(mxd, "", ryx):

if lyr.name == mapName:

arcpy.mapping.RemoveLayer(ryx, lyr)

print(mapName + "Has Done")

print("All Done")

5 结语

Python语言是当今最受程序设计语言之一，它是一种动态的、面向对象的脚本语言，是不受局限、跨平台的开源泉编程语言，功能强大，简洁易懂，应用广泛，不仅与Web、Internet联合开发和软件开发等，而且在科学计算、统计分析和教育方面发挥了独特的作用。本研究在Arcgis软件系统中应用Python语言解决了森林资源数据的年度更新，取得了满意的效果。