D O I：10. 3969/j.issn. 1672-2469. 2019. 01. 024
基于A rcG IS二次开发的
张青，胡成
(中国地质大学（武汉）环境学院，湖北武汉4300741)
摘要：文章以A r cG IS为平台，以P yth on语言为二次开发语言，利用A r cG IS强大的空间处理分析功能，绘制和处理研究区的多边形边界，并生成符合Tough2规则的泰森多边形网格，再根据数值模拟项目的需要提取网格信息写入到M E S H文件中，该文件可直接作为Tough2的输入文件，增加了网格文件生成的多样性和灵活性，提高了工作效率。
关键词：A rcG IS;Python;Tough2；泰森多边形
中图分类号：P641.2 文献标识码：A文章编号：1672-2469(2019)01-0087-05
地下热流及水流数值模拟程序1ugh2发布于1991年，经过20多年的开发，已发展为现今的新流体模块，
应用范围极其广阔，如C O：地质封存、地下水流模拟、环境影响评价、地热储藏工程等。虽然T u g h2可以模拟各类情况，但是推广上却较为困难；需要用户具备一定的编程基础[1-]。在模拟过程的初期，模型构建和网格剖分工作量巨大，后人开发出相应的网格生成工具，如MeshMaker、PetraSim、M V i e w等，但是这些网格生成工具并不具有地理信息的概念，也缺乏对网格图形强大的处理能力。
A r G I S是具有强大空间分析与处理功能的可视化软件，可以对具有空间信息的多边形进行编辑，其平台提供的二次开发P y t h o语言是一种跨平台的开源脚本语言，具有简单易学、开发效率高、伸缩性和可移植性高、可嵌人和稳定成熟等优势[4-]。自ArcGIS 9. 0引人Python作为脚本语言发展至今，Python已成为一种运用到A r c G I S中实现其数据管理与分析等多种功能自动化语言，Pytho 可以创建在独立脚本环境、工具箱或模型中运行的脚本。同时，Python可以轻松调用A r c G I S中各种工具，或访问其他人员开发的P y t h o模块[6-]，从 而实现程序化处理复杂数据，以减少工作量和提高工作效率。本文根据数值模拟项目实际需要，在 ArcGIS平台通过Python实现Tough2多边形网格的生成，为Tough2的进一步模拟提供便捷的基础。1方法与原理
1.1 Tough2网格文件规则
T o g h2模拟器需要的网格文件包括两部分，E L E M E单元文件和C O N N E连接文件，前者反映每个网格单元的基本信息，如单元标识、X Y Z坐标、参数分区标识等；后者反映每个单元网格之间的关 系，如交界面的面积、网格中心到交界面的距离、重力方向与网格中心连线的余弦值、热交换因子等，具体参数见表1和表2。其实，大部分参数都是网格固有的属性，但有些参数是需要根据实际情况进行赋值的，如参数分区标识和热交换因子。实 际上T〇ugh2支持任意形状的网格剖分，无论是三角性、矩形还是任意多边形都是可以的。
1.2泰森多边形生成原理
在一定范围内生成离散点，利用离散点构建三角形网，并作出每个三角形各边的垂直平分线，每 个离散点被各边垂线形成的一个多边形包围，这个多边形即为泰森多边形[8-]。泰森多边形网格非常符合T o g h2模拟器的特点，其不规则三角网的建立是根据Delaunay准则生成，在网格剖分上更具零活性，A r G I S平台中自带泰森多边形生成工具，
收稿日期：2018-08-01
作者筒介：张青（1993年一），男，硕士。
通信作者：胡成（1976年一），男，副教授。
糖果屋历险记2
表 1 E L E M E单元文件参数说明
1RO C KS3.0473E+060.08.4695E+047.6766E+04-1.000E+01单元标识材质名称体积热交换面积X坐标Y坐标Z坐标
表 2 C O N N E连接文件参数说明
8K HHJ—1—*—!—*—&—*—/—*—4—*—6—*—-—*—5
1213.7649E+023.7649E+026.0395E+030.0
相邻单元标识渗透性序号距离d1距离d2交界面积余弦值
利用P y t h u语言可轻松调用该工具，并对生成的多边形进行信息提取和图形分析，如图1所示。
图1 A r c G I S随机点生成泰森多边形
1.3网格文件信息获取
生成泰森多边形后可利用A r c G I S提供的内置 函数提取网格文件所需要的坐标、体积、面积等参 数，并可以利用P S h u基本函数计算距离、余弦 值等参数，再根据实际情况赋予属性参数，最后可 将提取和计算出的参数写入网格文件中。
2技术流程
E L E M E单元文件和C O H H E连接文件正是T u gh2所需的M E S H网格文件，先有E L E M E单 元文件，再有C O H H E连接文件；即先有网格单 元，再有网格关系。本次按先获取E L E M E文件 再获取C O H H E文件进行研究，技术流程如图2 所示。
图2技术流程图2.1 ELEME文件参数获取
(1%随机点的生成
在网格边界内生成一些随机点，通过随机点划 分出符合Delaunay准则的不规则三角网（T I H)，三 角形各边的垂直平分线即可形成泰森多边形的边。在A r c G I S中提供了 2种方法用于生成随机点，一 种是使用创建随机点的工具，可以按照工具给出的 提示生成指定范围内的随机点；另一^种是调用A p> GIS中CreateRandomPoints_ management")函数。在 该函数中只需要给定边界文件、随机点个数、最小 点间距等参数便可快速生成随机点，其函数的部分 参数说明见表3。
表3随机点函数参数说明
参数说明数据类型
out_path随机点文件输出位置Feature
out_name随机点文件名称String constraining_featnre_class边界范围文件Feature number—of_points—or_field随机点数量Long；Field minimum_allowed_distance最小点间距Field;Linear unit
生成的随机点文件中，系统已经默认赋予点文 件一些属性，如X编号、X坐标、6坐标等，但 这些系统分配好的参数并不能满足T u g h2网格文 件生成所需要的信息，还需在随机点文件中添加材 质名称这一属性。
(2)创建泰森多边形
创建泰森多边形的前提是在已知点位置的条件 下，根据Delaunay准则生成的。与创建随机点类 似，ArcGIS也提供了 2种方法用于创建泰森多边 形：一种是利用CreateThiessenPolygons工具箱；另一^种是用 CreateThiessenPolygons_ analysis()函数，同样该函数只需要输入随机点文件和输出位置便可 以生成泰森多边形，具体参数见表4。
表4泰森多边形函数参数说明
参数说明数据类型
in_features随机点文件Feature
out_feature_class输出位置Feature
与随机点函数不同的是，生成的泰森多边形不 需要指定输出范围参数，所以，生成的泰森多边形 是一个矩形，该矩形的范围是输人点输人要素的范 围另加10%，在生成后还需用8*工具剪裁至模 拟区边界形状。
(3 %泰森多边形网格信息提取
生成符合要求的泰森多边形后需要提取每个单 元格的体积，此次研究二维网格的生成，实际上单 元体积在数值上等于单元面积乘以模拟层厚度。X、Y、Z的坐标值则需要根据重力方向进行选取，默认重力方向在z坐标轴上，竖直向下，X坐标和 Y坐标等于投影坐标系减去投影带号所在范围的基 数值，Z坐标等于模拟层厚度的一半。如果重力方 向为Y方向，则Y坐标与Z坐标互换，Y坐标在 数值上等于模拟层厚度的一半。对于平行多层网 格，在空间上改变的只需是Z值。
在A p G IS中，单元格面积和坐标是图形固有 属性，上述网格信息可以通过数据访问模块(arcpy.da%中的Search Cursor类进行只读权限访问。
至此，ELEME文件的主要参数ID索引码、X 坐标、Y坐标、单元体积、参数获取完成。
2.2 CONNE文件参数获取
(1 %创建邻近关系表
A r G IS中的面领域工具用于创建领近关系表，面邻域工具按照等级路径确定要在输出表中记录的 邻域类型和统计数据。邻域关系、重合边和结点邻 域按照从高到低的等级顺序进行重叠。一旦到更 高顺序的邻域，该工具计算和存储关系信息并跳过 较低顺序关系的分析。由于在创建的泰森多边形网 格（T0UGH2网格）中，领域类型不存在重叠面，需要的是重合边以获取相交面积这一参数。
在图3中，左边输人要素是一个3 x3的网格 单元，单元编号从101 ~ 109，右边是通过面领域 工具分析后，输出的面领域分析表。
在图3输出的面领域关系表中，揭示了相邻网 格单元间的关系，与101网格相邻的有102、104、105网格单元，102和104网格单元的领域类型是 重合边，105网格单元的领域类型是结点领域。对
INPUT
POLYGONS
107108109
104105106
101102103
100m
OUTPUT TABLE
OBJECE)* src myCode nbr myCode LENGTH NODE COUNT
n:^101102100
2101104100
31011050
4102101100
5102103100
61021040
7102105100
81021060
9103102100
101031050
11103106100
12104101100
图3面领域关系表
于T〇u g h2网格文件需要获取的是重合边，在图5 中是LENGTH属性字段。输出领近关系表后通过 数据访问模块（arcpy.da)中的SeerchCursor访问这 些信息。
(2%获取公共边到节点的距离
每两个相邻单元格都只有两个共同的端点，通 过数据访问模块进行访问，比较每个网格端点的X 坐标和Y坐标是否相等来获取公共端点的X坐标 和Y坐标，这样每个单元格中心与两个端点间便构 成了一火
个三角形，可通过海伦公式：
S=sqrt(p(p-a)(p-h)(p-c['%
p= (a+7+c)/2
和面积公式S = 1/2a Z求出单元格中心到公共 边的距离9，即网格中心到公共边所形成等腰三角 形上的高Z。
mmpp
实际中泰森多边形具有相邻关系的两个单元格 中心到公共边的距离相等，通过两中心的坐标求出 距离，其距离的一半即公共边到节点的距离。
(3%获取 cosine值
在二维网格中如果取重力方向在Y坐标轴上，获取重力方向与连接线夹角的cosine值方法与计算 公共边到节点的距离方式类型，以相邻网格单元的 中心节点（随机点）构建直角三角形，cosine值等于 两中心节点间的高差比上两中心节点的距离。在构 建直角三角形前需要判断两种情况的出现，如果两 节点的X坐标相等，则c s in e值等于0，如果Y坐 标相等，则c s n e值等于1或-1。重力方向取X 方向与重力方向取Y坐标轴方向相似，在此不做说 明。如果重力方向在Z坐标轴上，则所有的cosine 值取 0。
至此，C0NNE文件的主要参数相邻单元格的 I索引码、单元节点距离、c s n e值的获取完成，其中，渗透性编号与网格单元的空间位置关系有 关，如果两中心节点在X平面上则为1，在Y平面 上则为2，在Z平面上则为3，其反映的是X Y Z方 向的渗透率，1、2、3只是用于识别渗透性方向的
代号，可以在网格文件生成时写入。
国家新型城镇化规划(2014-2020年)2.3 MESH文件参数整合
在获取Tough2中的ELEME和CONNE所需的 主要参数后，将这些参数按照指定的格式写入到 MEWH文件中，使用Python中的format函数格式化 输出，例如体积的格式化输出：Volunme = format (ELEME[i] [2]，'.4E f。使用 Python 中的 rjust 可以将字符进行右对齐，如：f qite(str(D1%.rjust (10))。为了将上述繁琐的步骤简单化，方便 Tough2网格文件的生成，以Python工具箱的形式 进行集成，其工具箱界面如图4和图5所示。
图 4 M E S H M a k e r-P o ly g o n,p yt 工具箱
图 5 F Ile-P o ly g o n•qo R 工具箱
图6的Python工具箱用于生成指定范围内的泰 森多边形网格（Tough2网格），图7的Python工具 箱用于输出MEWH文件。
24参数分区
由于在实际模拟项目中，研究区内各个网格单 元的参数是不同的，在进行模型识别时，不可能也 没有必要让每个单元都赋予不同的参数值，所以只 能进行参数分区，在每个参数分区内的各单元的参 数是相同的[10]。先在ArcGIS平台将研究区划分为 多个子区域，每个子区域内的单元网格参数一致，在ELEME单元文件中，每一行有一个位置记录该 单元的参数，用相应的代号表示，可以是英文字母也可以是数字或者组合，对于具有相同参数的单元 该代号相同。在ArcGIS平台中可利用Python语言 将这些落在同一个参数分区内的单元节点的参数字 段统一赋予相同的代号，并对同一参数分区内的网 格进行加密处理。
3方法应用
以安徽某一石油洞库模拟为例，为了预测在库 区施工及运营条件下，库区及周边地下水水位下降 情况及其对地质环境可能的影响，需要在分析研究 区水文地质条件的基础上，建立水文地质概念模 型。根据地质时代、成因、岩性及工程性质的不 同，将研究区的岩性分为三大类。其中研究区东南 部主要为第四系地层，成分以粘土、砂土、砾石为 主；中部主要为印支期花岗岩体，成分以钾长石、斜长石、石英、黑云母为主；西北部主要为寒武纪 片岩段，成分以云母石英片岩，顶部为千枚岩为 主；洞库建立在中部的花岗岩体区域。研究区参数 分区如图6所示。
南京信息工程大学学分制管理系统
在ArcGIS中勾出模拟区边界的矢量文件后，在Python脚本中输入需要生成的网格数以及网格中 心最小间距即可将模拟区剖分为泰森多边形网格，如图 7 所示 。
根据落入每个参数分区内节点赋予岩性代号，
并自动生成T o u g h 2需要的M E S H 网格文件。本次 处理的坐标系为“1954北京坐标系”，网格中心最 小间距为10:，网格数量为100个，网格层厚为 20：。网格文件结果如图5 —9所示。
T
R0CKS3. 0473E+06
R0CKS7. 0467E+06 R0CKS2. 1433E+07 R0CKS7. 4015E+06 ROCKS1. 1434E+07 R0CKS6.0084E+06 R0CKS7. 6194E+06 R0CKS4. 3228E+06 R0CKS3. 0491E+06 R0CKS5.6644E+068. 4695E+047. 6766E+04-1. 000E+01 8. 3942E+047. 6747E+04-1. 000E+01 8. 5438E+048. 3131E+04-1. 000E+01 8. 3552E+047. 8297E+04-1. 000E+01 8. 4287E+048. 2341E+04-1. 000E+01 8. 4231E+047. 9951E+04-1. 000E+01 8. 0549E+047. 9400E+04-1. 000E+01 8. 3611E+048. 1970E+04-1. 000E+01 8. 4773E+048. 0744E+04-1. 000E+01 8. 6464E+048. 1127E+04-1. 000E+01
图8生成的E L E M E 网格文件
编辑和分析功能。
(2)
本次编程是通过生成随机点再生成泰森多边
形，对于较为复杂的边界可以利用语言将局部 网格进行加密，提高模拟的精度。对于其他网格类 型，ArcGIS 也提供了渔网函数实现矩形网格剖分。
(3)
本次编写的P y t h u 工具箱，主要针对二维
多边形，并以泰森多边形网格的生成为主，可满足 众多数字模拟项目的需要。针对三维立体研究区， 需借助A r c G I S 的三维平台实现二次开发，有待进 一^步研究。
13. 7649E+023. 7649E+026. 0395B+03
0.015. 0387E+025. 0387E+027. 9324E+03 0.012. 7378E+022. 7378E+025.1227E+03 0.013. 4602E+013.4602B+017. 8487E+03 0.012. 1956E+022. 1956E+021. 8752E+04 0.013. 6323E+023. 6323E+027. 7631E+03 0.012. 0049E+022.0049E+021. 2885E+04 0.014. 6830E+024. 6830E+024. 9104E+O3
0.0
图9 生成的C O N N E 网格文件
在图8的E L E M E 单元文件中，第一列代表的 数值含义与表1中第一列的含义一致，即图5中的 数字1代表单元标识为1的网格单元；在图9的
C 0N N E 网格文件中，第一列第一行表示相邻网格
单元1和2的渗透性、距离等信息，分别与表2中 的说明—
对应。通过在A r c G I S 平台上生成
T o u g h 2所需的M E S H 网格文件，不仅提高了模拟
者的工作效率，也为今后T o u g h 2网格文件的生成 提供了新思路和新方法。
4结论
(1 )
在A r c G I S 平台中利用P y t h o n 工具箱可在
具有地理信息的实体数据上直接生成网格文件，可 对任意边界的研究区进行网格剖分，并实现网格的
C0NNE-
参考文献
[1]施小清，张可霓，吴吉春.T 0UGH 2软件的发展及应用[J ].工
程勘察，2009, 37(10): 29-34.[2 ] Pruvs K ，Oldenburg C M ，Moridis G J . T 0UGH 2 usePs guide ver sion 2[J ]. 1999.[3] 胡立堂，张可霓，邓媛媛.基于多边形的T 0UGH 2模拟器网格 生成方法[J ].工程勘察，2010, 38(12): 32-37.[4] 邵保华，田学志.谈Python 在Arcgis 地理处理中的应用[J ].林 业勘查设计，2012(02): 99-100.[5 ] Sanner M F . Pytlion ： a programming language for software integration
and development [ J ] . J Mol Graph Model ，1999，17 (01 ) ： 57-61.[6]朱道强，关海涛.Python 在A rcG IS 中的应用[J ].测绘与空间 地理信息，2013, 36(07): 152-154.[7 ] Kwon P . Spatial Data Process Using Pytiion S crijD t and Visual Pro gramming Language in ArcGIS [ J ] . Computer Technology and Appli cation ，2013，4(06).[8] 龙胤慧，廖梓龙.基于泰森多边形法的庆阳市面雨量计算[J ]. 河北工程大学学报（自然科学版），2012(03): 64-67.[9] Han D ，Bray M . Automated Thiessen polygon generation [ J ] . Water
resources research ，2006，42 ( 11 ).[10] 陈陆望，何建东，施小平，等.松散承压含水层水文地质参数 分区及水流场数值模拟[J ].现代地质，2015, 29 (04): 967-974.
质壁分离
(上接第86页）建成之后，2018年各个洪涝风险指 标均呈现一^定幅度的降低。
虽然城市建设引起的土地利用变化和局部地坪 高程变化会引起洪涝风险的变化甚至加剧，但是防 洪
排涝工程的建设运行加快了涝水入江入河、洪水 外排的速率，将洪涝风险消弭到最低。所以，尽管 人类建设活动造成了相对低地形成，但防洪排涝工 程发挥的巨大效用使得洪涝风险较大程度下降。因 此，“规划”和实施系统的防洪排涝工程对降低城 市洪涝风险具有重要意义。
参考文献
[J ].水科学进展，2016, 27(04): 485-491.[2] 程晓陶，李超超.城市洪涝风险的演变趋向、重要特征与应对方 略[J ].中国防汛抗旱，2015, 25(03): 6-9.[3] 国家防汛抗旱总指挥部，中华人民共和国水利部.中国水旱灾 害公报2016[M ].北京：中国地图出版社，2017(07).[4] 侯春生.城市洪涝灾害的形成及防治措施[J ].环境与发展， 2017, 29(03): 101-102.[5] 丁曼，李航，于得万.M IKE 模型在月亮泡蓄滞洪区洪水演进中 的应用[J ].水利规划与设计，2018(09): 49-51.[6] 刘艳红，靳志伟.二维水力学模型在红光大桥洪水影响评价中 的应用[J ].水利规划与设计，2017(11): 157-161.[7] 方国华，金德钢，刘俊，等.甬江流域河网水利计算模型研究及
[1]张建云，王银堂，贺瑞敏，等.中国城市洪涝问题及成因分析
应用[J ].灾害学，2009, 24(02): 41-45.
豫公网安备41160202000603
站长QQ:729038198
关于我们
投诉建议