鱼C论坛

 找回密码
 立即注册
查看: 3595|回复: 4

[技术交流] ARCGIS与python(1/3)

[复制链接]
发表于 2016-10-1 15:47:49 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能^_^

您需要 登录 才可以下载或查看,没有账号?立即注册

x
目录

前言
一、PYTHON  语言基础
1 数学运算符
2 字符串操作
3 模块的使用(M ODULES )
4 使用 DEF 构建函数
5 流程控制结构:I F ,W HILE ,F OR
6 简单输入和输出
二、ARCGIS&PYTHON
1 如何创建地理处理对象( GEOPROCESSOR OBJECT )
2 获取地理处理帮助
2.1 举例:如何使用 Geoprocessor Programming Model 中的 Lists
3 使用地理处理工具——T OOLBOXES 和 A LIASES
4 在建模中使用脚本(S CRIPTS IN  M ODEL B UILDER )
5 在 P YTHON W IN 里调试地理处理脚本
5.1  调试选择和消息
5.2PythonWin 的调试工具
5.3 地理处理工具举例
6 使用描述(D ESCRIBE )和存在(E XISTS )获取数据信息
6.1 描述
6.2 存在( Exists )
6.3 在循环中使用描述和存在
7 在 P YTHON 脚本中使用地图代数(M AP  A LGEBRA )
8 数据管理和指针(D ATA  M ANAGEMENT AND  C URSORS )
8.1 数据管理( Data Management )
8.2 指针( Cursors )
附录 1 :地理处理脚本中输入& 输出方法指南
附录 2 :其他

一直想学习 ArcGIS 中的 Python 脚本,大四下半学期终于有了时间,可是想找到这么一
本好的教材不容易。茫茫互联网,终于找到了旧金山州立大学 Jerry Davis 教授的个人主页,
对其中《Geoprocessing Scripts With Python》如获至宝,独乐乐不如众乐乐,现在将其教程
翻译并结合自己的学习情况给出总结。希望能够给更多想学习 Python 的同学一个参考。
另外,在我刚开始接触 Python 时,是看了台湾辅仁大学一位老师的视频课件,在此致
谢。
我想从两个大部分总结:一、Python 语言基础;二、ArcGIS&Python。其中第一部分参
考了《Python 精要参考(第二版) 》 、 《Python 编程金典(读书笔记) 》等书籍文献。对于多
数读者来说,可能或多或少有一些编程基础,所以理解起来应该不成问题。
文中多数数据来自 Jerry Davis 教授的主页,放在“C:\prog”目录下,为了直观,我将运
算结果一并编辑,方便参考。
值得一提的是 ArcGIS 的在线帮助文档, 一个实时更新的 GIS 宝库, 很多专业性知识都可
以找到答案,点击链接 ArcGIS10 中文帮助、ArcGIS9.3.1 或 9.3 英文帮助。 获取更过脚本例
子来学习 :ESRI 的地理处理模型和脚本工具库是个不错的选择。
由于我也是初次接触,翻译或者心得难免有纰漏之处,希望同仁们可以多多交流!
前言
在 GIS 建模或 GIS 数据管理中,你可能经常需要处理一系列步骤才可以完成的工作;你
可能有一个工作目录下的数据需要重投影、 裁剪到研究区域, 或者用某种方法组合成期望的
结果;我们也经常需要根据不同情形用不同方法处理数据,因此我们需要作出选择,而高质
量的决策需要考虑很多低水平的决策,这可以通过脚本程序模型辅助完成。
脚本编程的主要目的是使枯燥的处理数据工作自动化, 通过逻辑来指挥处理过程。 我想
自动化和逻辑是关键,它们区别于我们多数使用计算机时的交互活动。我们发 E-mail,写文
章或者设计地图,都需要和计算机交互,而处理一系列数据,我们需要自动化和利用逻辑来
指导自动化。
在地理处理脚本逻辑中,我们需要在允许我们做的事情中作出决定,比如,处理栅格数
据不同于矢量数据,或为没投影的数据设置投影,或处理仅在特定时间搜集的数据集。对于
重要的 GIS 工作来说,脚本以及其他形式的程序是必需的,而非可有可无。
在接下来的联系中,我们会探索 Python 的使用以及创建脚本来使用 ArcGIS 里众多的地
理处理工具。所有你能在 ArcToolbox 或 Model 中使用的工具都能够用在 Python 脚本中,这
些脚本可以生成脚本工具,像其他地理处理工具一样使用。

一、Python  语言基础
装 安装 PythonWin,在…\ArcGisDesktop9.3.iso\Desktop\PythonWin 目录下可以找到
PythonWin 的安装程序,默认是 不 安 装 的 ,
。同时会安装 win32com 以及
允许任何脚本在基于 Dispatch
的 地 理 处 理 过 程 中 工 作 。
ArcGIS10 中 引 入 了 全 新 的
Python Window 来增强内嵌的
Python 体验。
警告:不要尝试更新随
ArcGIS 安装的Python 到一个新
的版本!
下面介绍Python的一些简
单语法和规则。
1  数学运算符
Python 提供了多样化的通用数学运算符——多数编程语言的特征,以及许多通过
import 的 modules 提供的符号。常用的有+,-,*,/,**(幂),%(取模,即除后的余数)。
下面的表格显示了整型(Integer)和浮点型(Float)各种组合运算的结果,记住一
条规则,只要参与运算的有浮点型,则结果为浮点型;全为整型时,结果才为整型。
4 / 33
输入表达式  结果 Notes
2+3  5  整型结果
2.+3  5.0  2.是浮点型,结果浮点型
2-3  -1
2*3  6  整型结果
2.*3  6.0  浮点型
5/2  2  整型
5./2  2.5
5%2  1  取模
Az=270
Newaz=az+180
Print newaz%360
90  取模的用途之一——方位角加 180 后逆转方向
5**2  25
25**0.5  5.0  没有 sqrt()功能,除非添加 math 模块
2  字符串操作
用 注:使用 Python  帮助:有超过 30  种内置方法来处理字符,请到 Sequence Types
下的 String Methods  寻找帮助!
字符串是一串字母,比如’San Francisco’,字符串下标从 0 开始。学习字符串语法的
最好方法是自己动手尝试,下标展示之:
输入  结果  Notes
print 'zhulj'.capitalize()
Zhulj  s.capitalize()即将 capitalize()方法用于 s  s='zhulj'
print s.capitalize()
print s[0]  z
Strings 可以像一个字母列表一样处理,第一个字
母下标为 0, 某个字符段可以用 1:3 来格式化: 从
第 1 个的开头到第 3 个的开头,不包括下标为 3
的字母; s[-1]表示倒过来第一个, 相当于 s[len(s)-1]
s1=s[1]
print s1
h
print s[-2:]  lj
print s[2:3]  u
print s[2:4]  ul
print s[2:],s[:5]  ulj zhulj
s2=s.upper()
print s2
ZHULJ  我们可以将字符串方法的结果赋给新的变量
s3=s+s2
print s3
zhuljZHULJ  字符串组合用“+”
print s*3  zhuljzhuljzhulj  字符串重复用“*” ,后为重复次数
selstr='"elev">1000'
print selstr
"elev">1000
字符串可以使用单引号或双引号,跨行时用双引
号。  othersel=”’elev’>1000”
print othersel
‘elev’>1000
print s.isupper()  False  一些方法返回值为布尔型(True 或 False) ,一些
返回索引值(下标值)  print s2.isupper()  True
5 / 33
p='d:/work/lu.shp'
print p.find(‘.’)
10
print p.find(‘/’)  2
plist=p.split('/')
print plist
['d:', 'work',
'lu.shp']
你可以用 split()方法解析出不同的字符串片段, 并
创建一个列表(List) ,我们可以使用其中不同的
元素
print plist[0]  d:
print plist[1]  work
p2='d:\\work\\soil.shp'
print p2
d:\work\soil.shp
反斜线“\”和某些字母一起有特殊用法,如\n
为换行, “\”为转义字符,如“\\”则表示“\”
print 'Jerry\'s Kids'  Jerry's Kids
print 'Jerry\'s\nKids'
Jerry's
Kids
p3=r'd:\work\soil.shp'
print p3
d:\work\soil.shp
字符串前加“r”则强制“\”代表其本身,而非
转义字符,这对于文件路径的操作很方便
3  模块的使用(Modules )
Python 提供了一系列内置的方法(大量依赖于模块)用于通用编程。Python 安装时
自带了大量 Modules,最常用的有 math,sys,random,array 以及 os.path。
当然还有好多 Modules 可以下载,比如数字处理(Numeric)——numpy,可在
www.python.org 或 www.google.com 里搜索。 www.python.org/moin/NumericAndScientific
页面中列举了一些。
使用 Module 前,必须 import 之。通常我们会将一行 import <Module  名>放在程序
顶部,比如:
import arcgisscripting
当然,这不必成为你程序的第一行,但必须在使用它里面方法之前。
当要引用多个模块是,中间用逗号分隔,比如:
import arcgisscripting,sys,string,os,math
我们也可以自己为频繁使用的方法创建 Module,下面,我们开始体验内置的
Modules。
math 和 和 random  模块
很多常用的数学计算功能都可以通过 math 找到,比如三角计算或对数计算,如果
要使用复杂数字,就使用 cmath 模块。
和之前一样,通过以下表格来体现模块的使用:
输入  结果  Notes
import math
print math.log10(100)
2.0  以 10 为底的对数
print math.log(100)  4.60517018599  自然对数
print math.pi  3.14159265359  π 是一个静态常量,所以不需要括弧
6 / 33
pi=math.pi
print pi
3.14159265359
如果不想总是输入“math.pi”可以将
其赋给一变量
pi  3.1415926535897931  不需要 print 即可查看变量值
print math.sin(radians)
print math.cos(radians)
print math.tan(radians)
三角函数的计算是弧度制
degrad=pi/180
45*degrad
0.78539816339744828  度转化为弧度
sin=math.sin
sin(45*degrad)
sin(90*degrad)
0.70710678118654746
1.0
即使功能函数(像 sin)都可以赋给一
个变量
math.e  2.7182818284590451
math.hypot(3,4)  5.0  此方法是求三角形的斜边
x1=520382;y1=4152373
x2=520475;y2=4152963
不同赋值语句间用“; ”分隔
xr=x2-x1
yr=y2-y1
math.hypot(xr,yr)
math.sqrt(xr**2+yr**2)
(xr*xr+yr*yr)**0.5
597.28468923956189  不同的方式,相同的结果
import random
random.random()
0.27281588185756478
random()方法,每次结果都不同,值
域为[0.0,1.0)
rnd=random.random
rnd()
0.4456393835072503
mu=50
s=10
print
random.gauss(mu,s)
46.528202194
4  使用 def  构建函数
有点像 Module,但更简单,函数是一个自己定义功能,用在之后的代码中,并且
提供任何你想要使用的参数。这个函数从此可像变量那样在程序中使用,结合例子更容
易理解。接下来的代码定义了一个将度转换为弧度的简单函数,同时也定义了一个弧度
转换为度的函数,它们和 Excel 内置的函数类似。
import math
def radians(angdeg):
return angdeg*math.pi/180
def degrees(angrad):
return angrad*180/math.pi
print math.sin(radians(45))
print degrees(math.acos(0.5))
运行之,得到结果:0.707106781187 60.0
7 / 33
5  流程控制结构:If ,While ,For
任何脚本或编程语言的一个重要特征就是执行一系列不同情形语句的能力。
你想要创建一系列山影栅格来代表夏天、冬天和春秋分。山影(hillshade)工具需
要有太阳高度角和方位角作为输入参数。
重要日期  太阳倾角
夏至(6 月 21 日)  23.44
春秋分(3 月 21 日,9 月 21 日)  0
冬至(12 月 21 日)  -23.44
接下来是一段相当简单的代码,通过太阳倾角(太阳光线正午垂直照射的纬度)获
取太阳角和方位角以及纬度。输入两个参数:lat(研究区域的纬度,南半球为负)和
decl(太阳倾角) ,由此得到 sunangle 和 azimuth:
lat=30
decl=20
sunangle=90-lat+decl
azimuth=180
if sunangle>90:
sunangle=180-sunangle
azimuth=0
print sunangle,azimuth
上面的例子中 lat 和 decl 强制赋了值。
有三种流程控制操作:
if 仅在一个特定情形下才执行语句;
while 当一种情形存在下,持续执行语句
for 遍历一系列值
这些语法和 def 有些相似:初始语句后加顿号、需要执行的语句块有缩进。
这三个结构的一些重要的公共特征:
①if、while、for 语句均以冒号结尾,接下来是缩进的代码块,用于 if、while、for 定义
的情形。在脚本编写窗口,你会发现,你在一行末尾打上冒号后,下一行自动缩进,在接下
来的一行按下退格键取消缩进。
②如果你只需做一件事情,你可以在冒号后面同一行添加简短的语句,比如:
if x>0: print ‘x 比 0 大’
print ‘下一行不要缩进了。 。 。’
if (continued )
接下来,我们会探索一下另一个方便的模块:os.path:
开始之前,在 d:/下创建一个“testfolder”文件夹,然后新建一个“test.txt”文件;
尝试以下代码段,确保 print 语句前有缩进。
import os.path
if os.path.exists("d:/testfolder/test.txt"):
print "测试文件夹存在"
print "txt 文件存在"
8 / 33
elif os.path.exists("d:/testfolder"):
print "测试文件夹存在"
print "测试文件夹存在,但 txt 文件不存在"
else:
print "两者都不存在"
可选探索示例
接下来的例子做的事情对 GIS 非常重要,但是实际上不用任何地理处理代码。USGS7.5
米分辨率 DEM(数字高程模型)是文本文件(USGSDEM 文件) ,投影为 UTM,UTM 北向和
东向单位是米,但是高程单位可能是英尺(feet)或米(meters) 。因此在获取垂直或水平距
离信息时会有问题,比如坡度可以通过垂直距离/水平距离获得。如果你不在使用 Z 值之前
设置为 0.3048,将会出现错误结果。但是不幸的是,你可能不知道 DEM 文本文件的垂直单
位是英尺还是米。这些信息保存在第 539 个字符里, “1”代表英尺, “2”代表米,所以可以
通过读取这个文件判断。
下面的脚本演示了上述内容:
import fileinput
infile=r"c:\prog\pendata\woodside.dem"
firstline=fileinput.input(infile)[0]
unitchar=firstline[539]
unit="(unknown:not a 7.5' DEM?)"
if unitchar=="1":unit="feet"
if unitchar=="2":unit="meters"
print "\nElevation in"+" "+unit
fileinput.close()
输出结果:Elevation in feet
while (continued )
&#61656;  运行下面的代码,说明了一种 while 循环:
x=1
while x<10:
print x
x=x+1
屏幕依次输出 1~9
&#61656;  下面说明一下“==” (等于)的概念:
x=5
z=x==4
print z
输出 False
x=5
z=x==5
print z
输出 True
“==”是逻辑运算符之一,其他有“<” (小于) 、 “>” (大于) 、 “>=” (大于或等于) 、
“<=” (小于或等于) 、 “<>” (不等于) 。
使用逻辑运算符计算得到结果为布尔型:true(1)和 false(0) 。
下面例子简单体会一下布尔型表达式:
x=1
while x<10:
9 / 33
print x
x=x+1
表达式“x<10”结果是 true 或 false,所以这样允许我们在计算完一种情况时运行
一系列代码。许多情况下我们需要使用条件代码。
while 循环的一个优点是允许我们跳过整个部分,如果条件不满足初始情况。由于
while 提供一种容易结束循环的方法,我们甚至用它代替 if 语句。当循环一个数据
集时(GIS 中很常用的工作)while 循环很有用。在后面地理处理中我们会接触一些
例子。
for
&#61656;  尝试下面代码,演示了 for 循环:
for x in [1,2,3,4]:(注:[1,2,3,4]用 range(1,5)代替是一样的)
msg="hello world"
print str(x)+" "+msg(注:当我们希望一个数字 x 和字符串连接时,必须先
对数字进行格式转换即 str(x),否则系统报错)
&#61656;  下面的代码创建并输出指定文件夹内 shp 文件名列表(每个都以‘.shp’结尾)
import os
ws="c:/prog/hmbarea"
ilist=os.listdir(ws)#创建一个列表保存工作文件夹内的文件
fcs=[]#创建一个空列表,保存结尾为‘.shp’的文件
for i in ilist:
if i.endswith(".shp"):
fcs.append(i)
for fc in fcs:
print fc(输出结果如右图所示)
&#61656;  下面这个例子的循环较多次数。变量 mu 是算术平均数,s 是标准差——这两个是
random.gauss()用到的参数,可以尝试改变 n 的值查看结
果!
import random
mu=50
s=10
z5=mu-s*1.96
z95=mu+s*1.96
count=0
n=100000
for i in range(n):
x=random.gauss(mu,s)
if x<z5 or x>z95:count=count+1
print float(count)/n(每次运行的结果都不同,但都在 0.05 左右,即统计结果
在 5%左右)
6  简 单输入和输出
我们现在利用前面计算太阳角代码的片段,之前是直接指定参数值,现在我们有很多种
10 / 33
方法提供输入参数,现在我们用 sys 方法。尝试下面的代码,点击 (run 运行)之
后,在对话框中函数自变量(Arguments)一栏填入:40 23.44,如下图:
import sys
lat=float(sys.argv[1])
decl=float(sys.argv[2])#使用 arguments(argv)方法从“Arguments”一栏中获取输入
参数,并指定一个浮点型转换将字符型输入值传递给 lat 和 decl
sunangle=90-lat+decl
azimuth=180
if sunangle>90:
sunangle=180-sunangle
azimuth=0
print "正午太阳角="+str(sunangle)
print "方位角="+str(azimuth)(结果:正午太阳角=73.44 方位角=180)
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
回复

使用道具 举报

发表于 2016-10-2 12:50:12 | 显示全部楼层
转载要给源地址

代码要用 代码格式
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
回复 支持 反对

使用道具 举报

发表于 2017-2-4 19:06:06 | 显示全部楼层
谢谢分享!
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
回复 支持 反对

使用道具 举报

发表于 2019-5-10 11:44:54 | 显示全部楼层
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
回复

使用道具 举报

发表于 2022-6-27 17:21:57 | 显示全部楼层
1
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
回复

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

小黑屋|手机版|Archiver|鱼C工作室 ( 粤ICP备18085999号-1 | 粤公网安备 44051102000585号)

GMT+8, 2024-11-17 18:32

Powered by Discuz! X3.4

© 2001-2023 Discuz! Team.

快速回复 返回顶部 返回列表