科班程序员转gis开发之路(2019.12.19):python矢量切影像(shp切tif),并且逐个导出用于切割的要素对象(shp里的面,含属性)

一、开发背景介绍

	之前都是做供应链项目的,哪曾想到现在要转gis的业务了。
	目前也就是做一些基本的数据处理的工作。
	周围的大哥们都是用C#开发的,就我一个用python写脚本处理数据(毕竟调用方便)。
	由于我也是初次接触gis的业务,很多专业术语我也不懂,我会尽量以文件的名字取指代对应的专业术语。
	多多指教

二、矢量切影像的方法

1、基于arcpy,即调用arcgic相应系统工具接口。

	这里用到的函数接口都是arcmap软件里的系统工具提供的,用python的原因很简单,因为我不会用arcmap 0.0
	用python的真正原因是为了自动化,根据我写的逻辑,多次调用系统工具,不用手去一个个点了。
	或者,我在生成文件的时候,自动的把生成的文件放到我规定的某个自动生成的文件夹下...也是会方便很多。
	这份代码运行的环境是2.7,要注意哦是2.7
	话不多说上代码:
# -*-coding:utf-8 -*-

import arcpy
import os
import sys
import re

#过滤非法字符,我这里只过滤了‘-’,并把它替换成了‘_’
def validateTitle(title):
	# "[\/\\\:\*\?\"\<\>\|]"  # '/ \ : * ? " < > |'这个写法我是测试过的,是可以用的
	#因为我这里面有-,不能命名成shp文件,我也不知道为什么,反正不能命名,所以我要过滤掉‘-’
	new_title = re.sub('[-]', '_', title)
	return new_title
#这5行代码是用来修改py文件‘默认’编码方式的(注意修改的是默认编码)
#虽然我们第一行声明了编码格式是utf-8,但是这是python2(前面让大家一定注意版本,坑就在这里)
'''(补充:我刚刚查了为什么要设置编码格式为utf-8,因为如果不设置,那么中文注释就会报错。
但如果,你在编辑器里设置了文件编码格式是utf-8,这个声明写不写都可以,
所以为了统一不同编辑器下设置,俺们程序员统一会在第一行写这个声明。
欸!我想起来,好多人是在linux上开发python,是必须得在开头声明这个,不然写注释都报错)
'''
#python2的默认编码是ascII码,虽然我们声明了编码格式为utf-8
#但是他的默认编码格式是不会改的,不信你打印sys.getdefaultencoding()看看,是不是ascII 0.0
#所以当我们想要把一个字符串转成其他编码的时候,比如str="12ss接口"	str.encode('gbk')
#这样的代码能运行吗?
#当然不行。因为我们开始声明了文件编码格式是utf-8,在python2中字符串的编码方式就是utf-8
#而encode()函数是将字符串重新编码成gbk格式的,这里有一个细节隐藏了,
#编码之间的转换只能由unicode编码格式作为中间量
#比如utf-8转换成gbk,必须得utf-8->unicode,unicode->gbk
#我们成任何编码转换成unicode叫解码,unicode转换成gbk叫编码(这个‘编’是动词哦!)
#懂我意思吧。编码格式的直接转换是不可以的~,但是为了方便大家使用,encode()做好事不留名
#帮大家做了个解码的过程也就是decode('utf-8'),举个栗子:
#str.encode('gbk')等价于str.decode('utf-8').encode('gbk')
#utf-8是str的编码格式,decode按照utf-8的格式将str解码成unicode
#unicode在encode加密成gbk格式
#str.decode('utf-8').encode('gbk')是我们先指明str的编码格式进行的解码,再编码
#str.encoded('gbk')我们没有指明str的编码格式,所以它就通过sys.getdefaultencoding()获得到默认的编码格式
#也就是ascii,是不是很想锤死它。str的编码格式是和文件编码一样的,也就是我们开头声明的utf-8。
#我用utf-8上的锁,你怎么想用ascii开啊?不配套啊!0.5的笔芯和0.7的铅笔是不会有好结果的 0.0
#所以我们就用下面的5行代码,把默认编码给改了,斩草除根。
#那么问题又来了,改设置嘛,直接用下面的sys.setdefaultencoding('utf-8')不就可以了吗?
#写这么多干什么啊?原因很简单,因为你直接写sys.setdefaultencoding('utf-8')编译器会报错
#说,我没有这个方法啊!!!
#???!!!
#因为这个方法在我们import sys完成以后,就没了,我也不知道为什么没了,反正就是没了
#得重新加载这个库就有了,也就是reload(sys)
#但是reload以后,sys里的print函数又会没了0.0
#所以,我们在import sys 以后,先把里面一些,可能受影响的模块给存到变量里
#然后reload,重新加载sys库
#再把现在的模块用之前暂存的替换掉
#事实上如果你不用encode()的话用不着这样做,根本不影响,我就闲的慌,顺便理清楚各种编码和解码之间的关系,也是拾人牙慧,做个分享。
stdi, stdo, stde = sys.stdin, sys.stdout, sys.stderr
reload(sys)
sys.stdin, sys.stdout, sys.stderr = stdi, stdo, stde
sys.setdefaultencoding('utf-8')
print sys.getdefaultencoding()

shpFile = r"C:\123\456\789.shp".decode('utf-8')
outPath = r"D:\a\b\c"
inTifFile = r"E:\d\e\f.tif"

i = 0
with arcpy.da.SearchCursor(shpFile, ["SHAPE@", "NAME"]) as cursor:
	for row in cursor:
		#生成以name字段的值为名的文件夹
		if os.path.exists(os.path.join(outPath, row[1])) is False:
			os.makedirs(os.path.join(outPath, row[1]))
		ss = "正在切割出第{}块tif".format(i)
		#将ss解码成unicode,unicode可以显示中文 0.0
		print ss.decode('utf-8')
		#mask就是掩膜的意思,切割出掩膜位置形状的影像
		mask = row[0]
		extent = str(mask.extent.XMin) + " " + str(mask.extent.YMin) + " " + str(mask.extent.XMax) + " " + str(mask.extent.YMax)
		"""
		参数分别是,被切割的影像文件
		掩膜的最小矩形块坐标(想象一下怎么把一个矩形框正好套在一个掩膜形状上面)
		切割出来的一个影像文件(指定了绝对路径和文件名)
		noData的值,设置为0,或者-9999都可以,你想设什么设置什么
		这个参数是为了不按照掩膜的矩形块切割,前面不是指定了
		掩膜的最小矩形块嘛,可以按照矩形块切,也可以按照,指定的mask的形状去切。
		我这里指定的是按照mask的形状去切,None就是按照最小矩形去切
		"""
		arcpy.Clip_management(inTifFile, extent, os.path.join(outPath, row[1] + '\\' + validateTitle(row[1]) + '.tif', mask, "0", "ClippingGeometry"))
		#切影像的时候,你也可以顺便把切的矢量数据也给分别导出、
		"""
		参数介绍:
		被切割的shp文件,
		切割出来的一个面的文件名称(用shp单独存)
		选择条件,我这切出来的面,切的是哪一块(NAME属性 等于 多少多少的 面)
		"""
		arcpy.Select_analysis(shpFile, os.path.join(outPath, row[1] + '\\' + validateTitle(row[1]) + '.shp', "NAME = '{}'".format(row[1]))
	i+=1
del row
del cursor

三、小技巧

这里再说个python 的小技巧,因为python主要都是调用别人写好的接口,一般我们都是传参数就能用了。但是有很多的参数都是可选的,比如现在func(param1, param2, param3),1是必选,2,3是可选。但是我想传1,2。参数3不想传。那你写func(1,2)可能就会报错,说你参数3没传,目前这个我也没看透是什么机制,但是如果你指定好参数的名称,就不会报错了,如,func(param1 = 1, param2 = 2),这里的param1和param2是接口文档里面写的参数关键字,这个坑我可是趟了好几次才弄明白的,你们知道,找到一个接口,高高兴兴的去用,结果说你哪个参数没传,说好的是可选参数,你却逼着我硬要写。(因为接口不熟悉,有些配置参数不知道该不该给,所以我选择了让它选择默认的,所以我就不传了)
我是参考这里写的链接

四、小结

本来还打算写下用gdal库怎么去裁切的,但是遇到了问题,裁切出来的影像太大了,应该是用错了接口,等我再查查看接口,再写。
这里用的接口是属于arcmap软件里的系统工具的接口,用来裁切影像和合并影像再显示上面存在问题,和原图相比,看上去有小白点,不像原图,但数据是没有问题的,就是显示上有点问题,所以大佬们都用gdal去操作影像 0.0

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值