第8章代码

8.1  获得要素对象的属性和方法调用

#coding=utf8
import arcpy

import os
import sys
import math
import datetime
import numpy
from types import FunctionType

def getattrinfo(geo):
    if hasattr(geo,"area"):
        arcpy.AddMessage("{0}={1}".format("Area",geo.area))
    pnt=geo.centroid
    arcpy.AddMessage("centroid,x={0},Y={1}".format(pnt.X,pnt.Y))

    pnt=geo.firstPoint
    geoList.append(arcpy.PointGeometry(pnt))
    arcpy.AddMessage("firstPoint,x={0},Y={1}".format(pnt.X,pnt.Y))

    arcpy.AddMessage("hullRectangle={0}".format(geo.hullRectangle))
    arcpy.AddMessage("length={0}".format(geo.length))
    arcpy.AddMessage("isMultipart={0}".format(geo.isMultipart))
    arcpy.AddMessage("partCount={0}".format(geo.partCount))
    arcpy.AddMessage("pointCount={0}".format(geo.pointCount))
    arcpy.AddMessage("type={0}".format(geo.type))
    if geo.type=="polygon":#callable(geo.getArea): #方法是否存在
        arcpy.AddMessage("GEODESIC 面积={0}".format(geo.getArea("GEODESIC")))
    elif geo.type=="polyline":
        pnt=geo.positionAlongLine(0.5,True).firstPoint #线中点
        geoList.append(arcpy.PointGeometry(pnt))
        arcpy.AddMessage("line centerPoint,x={0},Y={1}".format(pnt.X,pnt.Y))


    #boundary=geo.boundary()
    #boundary=geo.buffer(50)
    boundary=geo.convexHull()
    geoboundary.append(boundary)


def main():
    fields=("OID@","SHAPE@")
    with arcpy.da.SearchCursor(inFeature, fields) as cursor:
        for row in cursor:
            FID=row[0]
            arcpy.AddMessage("================FID={}".format(FID))
            geo=row[1]
            getattrinfo(geo)

inFeature=arcpy.GetParameterAsText(0)
outFeature=arcpy.GetParameterAsText(1)
outFeature1=arcpy.GetParameterAsText(2)
geoList=[]
geoboundary=[]
main()
#arcpy.CopyFeatures_management(geoList, outFeature)
arcpy.Select_analysis(geoList, outFeature)
arcpy.CopyFeatures_management(geoboundary, outFeature1)

8.2  导出几何（点，线，面）要素坐标

#coding=utf8
import arcpy

import os
import sys
import math
import datetime

def text_save(filename, data):#filename为写入CSV文件的路径，data为要写入数据列表.
    file = open(filename,'a')
    for i in range(len(data)):
        s = str(data[i])+'
'

        file.write(s)
    file.close()
def getPoint():
    pList=[]
    pList.append("Point")
    rows = arcpy.da.SearchCursor(infc,("OID@","SHAPE@"))
    # Enter for loop for each feature/row
    #
    for row in rows:
        # Create the geometry object
        #

        FID=row[0]
        feat = row[1]
        pnt= feat.getPart(0)
        arcpy.AddMessage("FID="+str(FID)+",X="+str(round(pnt.X,xsnum))+", Y="+str(round(pnt.Y,xsnum)))
        pList.append("FID="+str(FID)+",X="+str(round(pnt.X,xsnum))+", Y="+str(round(pnt.Y,xsnum)))

    text_save(txtFile, pList)

def main():

    desc=desc = arcpy.Describe(infc)
    ptype = desc.shapeType
    if ptype=="Point":
        getPoint()
        return

    pList=[]
    pList.append(ptype)

    # Create search cursor
    #
    rows = arcpy.da.SearchCursor(infc,("OID@","SHAPE@"))

    # Enter for loop for each feature/row
    #
    for row in rows:
        FID=row[0]
        feat = row[1]
        arcpy.AddMessage("FID=%i" % FID)
        pList.append("FID=%i" % FID)

        partnum = 0

        # Step through each part of the feature
        #
        for part in feat:
            # Print the part number
            #

            arcpy.AddMessage("Part %i" % partnum)
            pList.append("Part %i" % partnum)
            # Step through each vertex in the feature
            #
            i=0
            for pnt in feat.getPart(partnum):
                if pnt:
                    # Print x,y coordinates of current point
                    #

                    arcpy.AddMessage("i="+str(i)+",X="+str(round(pnt.X,xsnum))+", Y="+str(round(pnt.Y,xsnum)))
                    pList.append("i="+str(i)+",X="+str(round(pnt.X,xsnum))+", Y="+str(round(pnt.Y,xsnum)))
                else:
                    # If pnt is None, this represents an interior ring
                    #
                    pList.append("Interior Ring:")
                    arcpy.AddMessage("Interior Ring:")
                i=i+1
            partnum += 1
    text_save(txtFile, pList)


infc=arcpy.GetParameterAsText(0)
xsnum=arcpy.GetParameter(1)
txtFile=arcpy.GetParameterAsText(2)

main()

8.3  导入坐标，生成多段线和带孔的面多部件面

#coding=utf8
import arcpy

import os
import sys
import math
import datetime
import codecs

def ReadTXT(txtFile):
    sumlist=[]
    f = open(txtFile,"r")
    lines = f.readlines()

    for line in lines:
        curline=line.replace('
', '') #删除


        sumlist.append(curline)
    f.close()
    return sumlist

def createPolyline(pList):
    num=len(pList)
    if num<3:
        return
    cursor = arcpy.da.InsertCursor(outFeature, ("SHAPE@"))
    array = arcpy.Array()
    sumarray= arcpy.Array()
    xystr=pList[1]
    FID=int(xystr.replace("FID=",""))
    for i in range(3,num):

        xystr=pList[i]

        if xystr.startswith("FID="):
            if sumarray.count==0:
                cursor.insertRow([arcpy.Polyline(array)])
            else:
                sumarray.add(array)
                cursor.insertRow([arcpy.Polyline(sumarray)])


            array.removeAll()
            sumarray.removeAll()
            FID=int(xystr.replace("FID=",""))
        elif xystr.startswith("Part "):
            #pass
            if array.count>0:
                sumarray.add(array)
                array.removeAll()
        else:
            xyList=xystr.split(",")
            n=len(xyList)
            if n>2:
                x=float(xyList[1].replace("X=","").replace(" ",""))
                y=float(xyList[2].replace("Y=","").replace(" ",""))
                #pnt = arcpy.Point(x, y,ID=FID)
                pnt = arcpy.Point(x, y) #加不加ID都可以
                #pnt=createonePoint(x,y)
                array.add(pnt)

    if array.count>0:
        if sumarray.count==0:
            cursor.insertRow([arcpy.Polyline(array)])
        else:
            sumarray.add(array)
            cursor.insertRow([arcpy.Polyline(sumarray)])


        array.removeAll()
        sumarray.removeAll()
    del cursor
def createPolygon(pList):
    num=len(pList)
    if num<3:
        return
    cursor = arcpy.da.InsertCursor(outFeature, ("SHAPE@"))
    array = arcpy.Array()
    sumarray= arcpy.Array()
    xystr=pList[1]
    FID=int(xystr.replace("FID=",""))
    for i in range(3,num):

        xystr=pList[i]

        if xystr.startswith("FID="):
            if sumarray.count==0:
                cursor.insertRow([arcpy.Polygon(array)])
            else:
                sumarray.add(array)
                cursor.insertRow([arcpy.Polygon(sumarray)])


            array.removeAll()
            sumarray.removeAll()
            FID=int(xystr.replace("FID=",""))
        elif xystr.startswith("Part "):
            #pass
            if array.count>0:
                sumarray.add(array)
                array.removeAll()
        else:
            xyList=xystr.split(",")
            n=len(xyList)
            if n>2:
                x=float(xyList[1].replace("X=","").replace(" ",""))
                y=float(xyList[2].replace("Y=","").replace(" ",""))
                #pnt = arcpy.Point(x, y,ID=FID)
                pnt = arcpy.Point(x, y) #加不加ID都可以
                #pnt=createonePoint(x,y)
                array.add(pnt)

    if array.count>0:
        if sumarray.count==0:
            cursor.insertRow([arcpy.Polygon(array)])
        else:
            sumarray.add(array)
            cursor.insertRow([arcpy.Polygon(sumarray)])


        array.removeAll()
        sumarray.removeAll()
    del cursor
#有坐标生成点
def createxyPoint(x,y):
    point = arcpy.Point() #create an empty Point object
    point.X = x
    point.Y = y
    #return point
    pointGeometry = arcpy.PointGeometry(point)
    return pointGeometry
#有坐标生成点
def createonePoint(x,y):
    point = arcpy.Point() #create an empty Point object
    point.X = x
    point.Y = y

    return point
def createPoint(pList):
    num=len(pList)
    cursor = arcpy.da.InsertCursor(outFeature, ("SHAPE@"))
    for i in range(1,num):
        xystr=pList[i]
        xyList=xystr.split(",")
        x=float(xyList[1].replace("X=","").replace(" ",""))
        y=float(xyList[2].replace("Y=","").replace(" ",""))
        pnt=createxyPoint(x,y)
        cursor.insertRow((pnt,))
    del cursor

def main():
    pList=ReadTXT(txtFile)
    num=len(pList)
    if num<1:
        return
    pType=pList[0].upper() #"后面后回车
    outPath, outFC = os.path.split(outFeature)
    arcpy.AddMessage("pType="+pType+",outPath="+outPath+",outFC="+outFC)
    arcpy.CreateFeatureclass_management(outPath, outFC, pType)

    if pType=="POINT":
        createPoint(pList)
    elif pType=="Polygon".upper():
        createPolygon(pList)
    else:
        createPolyline(pList)
txtFile=arcpy.GetParameterAsText(0)
outFeature=arcpy.GetParameterAsText(1)

main()

8.4.1       使用numpy导出点坐标

#coding=utf8
import arcpy

import os
import sys
import math
import datetime
import numpy
def text_save(filename, data):#filename为写入CSV文件的路径，data为要写入数据列表.
    file = open(filename,'a')
    for i in range(len(data)):
        s = str(data[i])+'
'

        file.write(s)
    file.close()


def main():
    arr = arcpy.da.FeatureClassToNumPyArray(infc, ("OID@", "SHAPE@XY"))
    num=arr.size
    pList=[]
    ss=str(xsnum)
    for i in range(num):
        yy="{0},{1:."+ss+"f},{2:."+ss+"f}"
        arcpy.AddMessage(yy)
        mystr=yy.format(arr[i][0],arr[i][1][0],arr[i][1][1])
        pList.append(mystr)
    text_save(txtFile, pList)
infc=arcpy.GetParameterAsText(0)
xsnum=arcpy.GetParameter(1)
txtFile=arcpy.GetParameterAsText(2)

main()

8.4.2使用numpy导入文本生成点

#coding=utf8
import arcpy

import os
import sys
import math
import datetime
import numpy
def ReadTXT(txtFile):
    sumlist=[]
    f = open(txtFile,"r")
    lines = f.readlines()

    for line in lines:
        curline=line.replace('
', '') #删除


        sumlist.append(curline)
    f.close()
    return sumlist
def main():
    sumlist=ReadTXT(txtFile)

    num=len(sumlist)
    idList=[]
    xList=[]
    yList=[]

    for i in range(num):
        mystr=sumlist[i]
        arr=mystr.split(",")
        #arcpy.AddMessage("==="+arr[0]+";"+arr[1]+";"+arr[2])
        idList.append(int(arr[0]))
        xList.append(float(arr[1]))
        yList.append(float(arr[2]))
    numpy_arr = numpy.rec.fromarrays([idList, xList, yList],
                              names=["编号","x","y"])
    if sr!=None:
        arcpy.da.NumPyArrayToFeatureClass(numpy_arr, outFeature, ['x', 'y'],sr)
    else:
        arcpy.da.NumPyArrayToFeatureClass(numpy_arr, outFeature, ['x', 'y'])

txtFile=arcpy.GetParameterAsText(0)
outFeature=arcpy.GetParameterAsText(1)
SR=arcpy.GetParameter(2) #坐标系
if str(SR)=="":
    sr = None
    
else:
    sr = arcpy.SpatialReference()
    sr.loadFromString(SR)
    
main()

8.5  线按指定长度分割

#coding=utf8
import arcpy
import sys
import os
import math
##########
def getCount(inFeature):
    result = arcpy.GetCount_management(inFeature)
    count= int(result.getOutput(0))
    return count

def initProgress(hint,num):
    arcpy.SetProgressor("step", u"正在"+hint,0,num,1)
def step():
    arcpy.SetProgressorLabel(u"正在进行....")
    arcpy.SetProgressorPosition()
def freeProgress():
    arcpy.ResetProgressor()


def splitNgeometry(sumgeometry,DIS):
    linedis=sumgeometry.length
    if linedis<DIS:
        return
    num=int(linedis/dis)
    for i in range(num):
        pnt=sumgeometry.positionAlongLine(DIS*(i+1), False)
        pointList.append(pnt)

arcpy.env.overwriteOutput = True

inFeature = arcpy.GetParameterAsText(0) #输入要素
dis = arcpy.GetParameter(1) #距离
outFeature = arcpy.GetParameterAsText(2)

jfb_Select=arcpy.env.scratchWorkspace+"/yl_999888"#不能c:要c:\或则 c:/
arcpy.MultipartToSinglepart_management(inFeature,jfb_Select)

pointList=[]
count=getCount(jfb_Select)
initProgress(u"处理",count)
try:
    rows = arcpy.da.SearchCursor(jfb_Select,"SHAPE@")
    for row in rows:
        step()
        geometry = row[0]
        splitNgeometry(geometry,dis)


    freeProgress()
    num=len(pointList)
    if num<1:
        arcpy.Select_analysis(jfb_Select, outFeature)
    else:
        arcpy.SplitLineAtPoint_management(inFeature,pointList, outFeature, "0.01 Meters")

    if arcpy.Exists(jfb_Select):
        arcpy.Delete_management(jfb_Select)

except Exception as e:
    arcpy.AddError(e.message)

8.6 要素类相对XY平移

#coding=utf8
import arcpy

import os
import sys
import math

def initProgress(hint,num):
    arcpy.SetProgressor("step", u"正在"+hint,0,num,1)
def step():
    arcpy.SetProgressorLabel(u"正在进行....")
    arcpy.SetProgressorPosition()
def freeProgress():
    arcpy.ResetProgressor()

##分解成多个列表，内部也分解
def splitNgeometry(mgeometry):
    num=mgeometry.count
    Sumarray = arcpy.Array()
    parray = arcpy.Array()
    for i in range(num):
        pt=mgeometry[i]
        if pt:
            parray.add(pt)
        else:#内边形
            Sumarray.add(parray)
            parray.removeAll()
    Sumarray.add(parray)
    return Sumarray

def MoveonePoint(pnt):
    x=pnt.X+moveX
    y=pnt.Y+moveY
    pt1 =  arcpy.Point(x,y)
    return pt1
def MovePoint(pt):
    pnt= pt.getPart(0)
    return MoveonePoint(pnt)
def MoveMPoint(partgeometry):
    array = arcpy.Array()
    num=partgeometry.partCount

    for i in range(num):

        pt=partgeometry.getPart(i)
        pt=MoveonePoint(pt)
        array.add(pt)
    return array

def MoveManyPoint(partgeometry):
    array = arcpy.Array()
    num=partgeometry.count

    for i in range(num):

        pt=partgeometry[i]
        pt=MoveonePoint(pt)
        array.add(pt)
    return array
def Moveobj(geometry):
    part_count = geometry.partCount #有几部分
    #arcpy.AddMessage("FID:"+str(FID)+",part_count:"+str(part_count))
    Sumarray = arcpy.Array()
    for i in range(part_count):
        partgeometry=geometry.getPart(i)
        SpliArray=splitNgeometry(partgeometry)
        N=SpliArray.count
        #arcpy.AddMessage("NNNNN=====:"+str(N))
        for j in range(N):
            Splitgeometry=SpliArray[j]

            array=MoveManyPoint(Splitgeometry)

            try:
                Sumarray.add(array)
            except Exception as err:
                arcpy.AddError("错误=============j:"+str(j)+","+err.message)
    return Sumarray
def Main():
    count=getCount(inFeature)
    if count < 1:
        arcpy.AddMessage(inFeature+u"没有数据")
        return

    arcpy.Select_analysis(inFeature, outFeature, '')
    if moveX==0 and moveY==0:
        arcpy.AddMessage(u"XY都为0，不用移动")
        return

    desc = arcpy.Describe(inFeature)
    layername=desc.name
    shapeType=desc.shapeType

    shapeName = desc.ShapeFieldName
    OIDField=desc.OIDFieldName
    initProgress(u"正在修改",count)


    rows = arcpy.UpdateCursor(outFeature)

    try:
        for row in rows:
            step()
            FID=row.getValue(OIDField)
            geometry = row.getValue(shapeName)
            #arcpy.AddMessage(u"================="+shapeType)
            if shapeType=="Point":
                geometry=MovePoint(geometry)
            elif shapeType=="Multipoint":
                #arcpy.AddMessage(u"多点=================")
                geometry=MoveMPoint(geometry)

            else:
                geometry=Moveobj(geometry)


            row.setValue(shapeName, geometry)
            rows.updateRow(row)
            arcpy.AddMessage("FID:"+str(FID)+"修改")
    finally:
        freeProgress()
        #del row
        del rows


##########
def getCount(inFeature):
    result = arcpy.GetCount_management(inFeature)
    count= int(result.getOutput(0))
    return count
inFeature  = arcpy.GetParameterAsText(0)
moveX  = arcpy.GetParameter(1)
moveY  = arcpy.GetParameter(2)
outFeature  = arcpy.GetParameterAsText(3)
Main()

8.7图形按属性合并（实现融合工具）

#coding=utf8
import arcpy

import os
import sys
import math
import datetime
import numpy


#获得字段类型是否为文本
def getFieldTypeIsText(TableName,FieldName):
    desc = arcpy.Describe(TableName)
    FieldName=FieldName.upper()
    for field in desc.fields:
        if field.Name.upper() ==FieldName:
            return field.type=="String"

    return False
def union(myvalue,isText):
    if isText:
        sql=inField+"='"+myvalue+"'"
    else:
        sql=inField+"="+str(myvalue)
    sumgeo=None

    with arcpy.da.SearchCursor(inFeature, ("SHAPE@"),sql) as cursor:

        for row in cursor:
            geometry=row[0]
            if sumgeo==None:
                sumgeo= geometry
            else:
                sumgeo=sumgeo.union(geometry)

    return sumgeo

def main():
    isText=getFieldTypeIsText(inFeature, inField)
    features = []
    with arcpy.da.SearchCursor(inFeature, (inField),sql_clause=(None,"group by "+inField)) as cursor:
        for row in cursor:
            myvalue=row[0]
            mygeo= union(myvalue,isText)
            features.append(mygeo)
    arcpy.CopyFeatures_management(features, outFeature)

inFeature=arcpy.GetParameterAsText(0)
inField=arcpy.GetParameterAsText(1)
outFeature=arcpy.GetParameterAsText(2)
main()

8.8  点按属性相同生成面

import arcpy
import os
import types

def convertPoints():
    arcpy.env.overwriteOutput = True

    # Input point FC
    # Output FC
    # Feature Field
    # Sort Field
    # Close Line or Leave Open
    inPts       = arcpy.GetParameterAsText(0)
    outFeatures = arcpy.GetParameterAsText(1)
    IDField     = arcpy.GetParameterAsText(2)
    sortField   = arcpy.GetParameterAsText(3)
    #closeLine   = arcpy.GetParameterAsText(4)

    if IDField in ["", "#"]: IDField = None

    if sortField in ["", "#"]:
        cursorSort = IDField
    else:
        if IDField:
            cursorSort = IDField + ";" + sortField
        else:
            cursorSort = sortField


    close = True

    convertPointsToLine(inPts, outFeatures, IDField, cursorSort, close)

def getZM(propType, hasMZ):
    envValue = getattr(arcpy.env, propType).upper()

    if envValue in ['ENABLED', 'DISABLED']:
        return envValue
    else:
        if hasMZ:
            return "ENABLED"
        else:
            return "DISABLED"

def convertPointsToLine(inPts, outFeatures, IDField, cursorSort, close):
    try:
        # Assign empty values to cursor and row objects
        iCur, sRow, feat = None, None, None

        desc = arcpy.Describe(inPts)
        shapeName = desc.shapeFieldName

        # Create the output feature class
        #
        outPath, outFC = os.path.split(outFeatures)
        arcpy.CreateFeatureclass_management(outPath, outFC, "POLYGON", "",
                                            getZM("outputMFlag", desc.hasM),
                                            getZM("outputZFlag", desc.hasZ),
                                            inPts)

        # If there is an IDField, add the equivalent to the output
        #
        if IDField:
            f = arcpy.ListFields(inPts, IDField)[0]
            fName = arcpy.ValidateFieldName(f.name, outPath)
            arcpy.AddField_management(outFeatures, fName, f.type, f.precision, f.scale, f.length,
                                      f.aliasName, f.isNullable, f.required, f.domain)

        # Open an insert cursor for the new feature class
        #
        iCur = arcpy.InsertCursor(outFeatures)

        # Create an array needed to create features
        #
        array = arcpy.Array()

        # Initialize a variable for keeping track of a feature's ID.
        #
        ID = -1
        fields = shapeName
        if cursorSort:
            fields += ";" + cursorSort

        for sRow in arcpy.gp.SearchCursor(inPts, "", None, fields, cursorSort, arcpy.env.extent):
            pt = sRow.getValue(shapeName).getPart(0)
            if IDField:
                currentValue = sRow.getValue(IDField)
            else:
                currentValue = None

            if ID == -1:
                ID = currentValue

            if ID <> currentValue:
                if array.count >= 2:

                    # To close, add first point to the end
                    #
                    if close:
                        array.add(array.getObject(0))

                    feat = iCur.newRow()
                    if IDField:
                        if ID: #in case the value is None/Null
                            feat.setValue(IDField, ID)
                    feat.setValue(shapeName, array)
                    iCur.insertRow(feat)
                else:
                    arcpy.AddIDMessage("WARNING", 1059, unicode(ID))

                array.removeAll()

            array.add(pt)
            ID = currentValue

        # Add the last feature
        #
        if array.count > 1:
            # To close, add first point to the end
            #
            if close:
                array.add(array.getObject(0))

            feat = iCur.newRow()
            if IDField:
                if ID: #in case the value is None/Null
                    feat.setValue(IDField, currentValue)
            feat.setValue(shapeName, array)
            iCur.insertRow(feat)
        else:
            arcpy.AddIDMessage("WARNING", 1059, unicode(ID))
        array.removeAll()

    except Exception as err:
        import traceback
        arcpy.AddError(
            traceback.format_exception_only(type(err), err)[0].rstrip())

    finally:
        if iCur:
            del iCur
        if sRow:
            del sRow
        if feat:
            del feat

        try:
            # Update the spatial index(es)
            #
            r = arcpy.CalculateDefaultGridIndex_management(outFeatures)
            arcpy.AddSpatialIndex_management(outFeatures, r.getOutput(0), r.getOutput(1), r.getOutput(2))
        except:
            pass

if __name__ == '__main__':
    convertPoints()