python编程结构（编写ArcGIS中的工具）--以我的第一个程序为例


"""
Tool Name:  Cost Path
Source Name: CostPath.py
Version: ArcGIS 10.2
Author: geying

This tool analyzes spatial patterns at multiple distances to assess whether the                       第一部分
observed pattern is more or less clustered/dispersed than might be expected
given random chance.  If an attribute field is specified, the field is                                            
说明：介绍整个程序的功能，主要包括红色字体部分
assumed to represent incident counts (a weight).

References:

Boots, Barry and Getis, Arthur. 1988. _Point Pattern Analysis_, Sage University

Paper Series on Quantitative Applications in the Social Sciences, series no. 07-001,

Beverly Hills: Sage Publications.

Getis, Arthur. 1984. "Interaction Modeling Using Second-order Analysis,

_Environment and Planning A_ 16:173-183.

"""

################### Imports ########################                                               ——     ###        ###    为注释

import os as OS

import sys as SYS

import collections as COLL

import numpy as NUM

import math as MATH                 

import arcpy as ARCPY                                                                                                                   第二部分

import arcpy.management as DM

import arcpy.da as DA                                                                                                                     导入模块

import SSUtilities as UTILS

import locale as LOCALE

LOCALE.setlocale(LOCALE.LC_ALL, '')

NODATA = 9999.0

RADIAN = 57.29578

KNIGHT_DIRS = [

    (-1, -1),                                                                                                                             ——  如果有常量，应写在这个位置

    (-1,  0),

    (-1,  1),

    ( 0, -1),

    ( 0,  1),

    ( 1, -1),

    ( 1,  0),

    ( 1,  1),

    (-2, -1),

    (-2,  1),

    (-1, -2),

    (-1,  2),

    ( 1, -2),

    ( 1,  2),

    ( 2, -1),

    ( 2,  1)

    ]

################ Output Field Names #################

outputFieldNames = ["NeighborhoodSize", "HeightMean"]

################### GUI Interface ###################

def setupCostPath():

    inRaster = ARCPY.GetParameterAsText(0)

    outputTable = ARCPY.GetParameterAsText(1)                                                                                      第四部分：创建输入输出

    neighborMethod = ARCPY.GetParameterAsText(2).upper()

    begSize = ARCPY.GetParameter(3)

    nIncrements = ARCPY.GetParameter(4)                                                          1、ARCPY.GetParameterAsText()  以文本的形式获取参数

    dIncrement = UTILS.getNumericParameter(5)

    saveIt = ARCPY.GetParameter(6)                                                                          ARCPY.GetParameterAsText().upper() 以文本的形式获取参数后将其全部转化成大写

    slopeDegreeLimit = UTILS.getNumericParameter(7)                                             ARCPY.GerParameter()  获 取参数                                                                                         

                                                                                                                                    UTILS.getNumericParameter() 获取double类型的参数

 #### Resolve Table Extension ####

    if ".dbf" not in OS.path.basename(outputTable):

        dirInfo = ARCPY.Describe(OS.path.dirname(outputTable))

        if dirInfo == "FileSystem":

            outputTable = outputTable + ".dbf"

    if neighborMethod == "ROOK":

        neighborMethod = "Rook"

    elif neighborMethod == "QUEEN":

        neighborMethod = "Queen"

    else:

        neighborMethod = "Knight"

    cp = CostPath(inRaster, outputTable, neighborMethod,                                               2、调用类     变量名 = CostPath(    ,     ,      ,     ,)   括号内的参数                                     

                  begSize = begSize, nIncrements = nIncrements,                                           与前面创建的输入输出相关联，通常可选参数用这种形式begSize = begSize表示将

                  dIncrement = dIncrement, slopeDegreeLimit = slopeDegreeLimit)  前面输入的值赋给一个名为begSize的变量（只是用了同样的名字而已，实际上具有不同的含义）

    cp.createOutput(outputTable, saveIt = saveIt)                                                            3、调用创建输出这个函数其规则为：变量名.函数名（参数1， 参数2）

class CostPath(object):

    """   """                                                                                                                                              第五部分：创建类

    def __init__(self, inRaster, outputTable, neighborMethod,

                 begSize = 3, nIncrements = 1, dIncrement = None, slopeDegreeLimit = 35):                       1、第一个函数通常为def__init__(self,      ,       ,      ,       ):

        #### Set Initial Attributes ####                                                                                                             后面的参数与前面创建的输入输出相关联，通常可选参数会给出缺省值

        UTILS.assignClassAttr(self, locals())

        ### neighboring process ###

        self.rook = self.neighborMethod == "Rook"

        self.queen = self.neighborMethod == "Queen"

        self.knight = self.neighborMethod == "Knight"

        #### Initialize Data ####

        self.initialize()                                                                                                                                 ——   函数的调用，程序运行到这行即跳转到类里的函数

        

        ###  Using relief function   ###

##        self.relief()

        ###  ###

        self.slope()

    def initialize(self):                                                                                                                                2、 def initialize(self):    是初始化函数，通常为类所需数值的初始化

        ### Create a Describe object from the raster band ###

        desc = ARCPY.Describe(self.inRaster)

        self.nrows = desc.height

        self.ncols = desc.width

        self.xSize = desc.meanCellWidth

        self.ySize = desc.meanCellHeight

        self.XMin = desc.extent.XMin

        self.YMin = desc.extent.YMin

        self.XMax = desc.extent.XMax

        self.YMax = desc.extent.YMax

        lower_left_corner = ARCPY.Point(self.XMin, self.YMin)

        self.lower_left_corner = lower_left_corner

        self.dem_mat = arcpy.RasterToNumPyArray(self.inRaster,

                                                lower_left_corner,

                                                self.nrows, self.ncols,

                                                nodata_to_value = NODATA)

       

        #### Determine All Distance Cutoffs ####

        cutoffs = []

        for inc in xrange(self.nIncrements):

            val = (inc * self.dIncrement) + self.begSize

            cutoffs.append(val)

           

        self.cutoffs = NUM.array(cutoffs)

       

    def relief(self):                                                                                                            3、定义了relief函数

        """     """                                                                                                                        函数下面通常有一段注释，用来说明函数的功能,用引号括起来

        

        #### Attribute Shortcuts ####                                                                                     空一行，整个函数也可以分为几个部分，用#号加注释

        M = self.nrows

        N = self.ncols

        Z = self.dem_mat

        M = 10

        N = 10

        ###  Create a zero matrix  ###

        myArray = NUM.zeros_like(Z)                                                                                  

       

        self.ld = COLL.defaultdict(float)

        ###  get sub_matrix calculate it mean ###

        for inc in xrange(self.nIncrements):

            BS = self.cutoffs[inc]

            k = BS / 2

            for i in xrange(M):

                iT = max(i-k, 0)

                iB = min(i+k, M)

                for j in xrange(N):

                    jL = max(j-k, 0)

                    jR = min(j+k, N)

                    sub_mat = Z[iT:iB+1 ,jL:jR+1]

                    hMin = NUM.min(sub_mat)

                    hMax = NUM.max(sub_mat)

                    myArray[i,j] = hMax-hMin

               

            hMean = NUM.mean(myArray)

            self.ld[inc] = hMean

        self.myArray = myArray

   

    def slope(self):                                                                                                             —— 定义了一个slope函数

        """ calculate the anisotropic accumualted costs """

        #### Attribute Shortcuts ####

        xSize = self.xSize

        Z = self.dem_mat

        M = 1

        N = 1

        slopeDegreeLimit = self.slopeDegreeLimit

        slope = slopeDegreeLimit / RADIAN

        maxSlopeDegree = 0

       

        self.slopeMat = NUM.zeros_like(Z)

        self.slopeDirMat = NUM.zeros_like(Z)

       

        fnSlopeRook   = lambda u, hp, ho: MATH.atan2(hp - ho, u)

        fnSlopeQueen  = lambda u, hp, ho: MATH.atan2(hp - ho, MATH.sqrt(2.0)*u)

        fnSlopeKnight = lambda u, hp, ho: MATH.atan2(hp - ho, MATH.sqrt(5.0)*u)

                   

        for i in xrange(M):

            for j in xrange(N):

                Ho = Z[i,j]

                ARCPY.AddMessage("Ho = %s" % Ho)

                for origID, origVal in enumerate(KNIGHT_DIRS):

                    ni, nj = (i + origVal[0], j + origVal[1])

                    ni = min(max(ni, 0), M )

                    nj = min(max(nj, 0), N )

                    Hp = Z[ni, nj]

                    if self.rook or self.queen or self.knight:

                        if origID in [1, 3, 4, 6]:

                            slope = fnSlopeRook(xSize, Hp, Ho)

                       

                    if self.queen == 1 or self.knight == 1:

                        if origID in [0, 2, 5, 7]:

                            slope = fnSlopeQueen(xSize, Hp, Ho)

                      

                    if self.knight == 1:

                        if origID in [8, 9, 10, 11, 12, 13, 14, 15]:

                            slope = fnSlopeKnight(xSize, Hp, Ho)

                    ARCPY.AddMessage("slope = %s " % slope)                                                              —— arcpy中的格式化输出

                    slope_degrees = slope * RADIAN

##                    if abs(slope_degrees) < slopeDegreeLimit:

##                        maxSlopeDegree = max(abs(slope_degrees), maxSlopeDegree)

##                       

##                self.slopeMat[i, j] = maxSlopeDegree

                       

                                       

                       

       

       

    def createOutput(self, outputTable, saveIt = False):                                                                 ——创建输出函数

        """ Creates Neighborhood Size Table  """

      

        #### Allow Overwrite Output ####

        ARCPY.env.overwriteOutput = 1

       

        #### Get Output Table Name With Extension if Appropriate ####

        outputTable, dbf = UTILS.returnTableName(outputTable)

       

        #### Set Progressor ####

        ARCPY.SetProgressor("default", ARCPY.GetIDMessage(84008))

        #### Delete Table If Exists ####

        UTILS.passiveDelete(outputTable)

       

        #### Create Table ####

        outPath, outName = OS.path.split(outputTable)

        try:

            DM.CreateTable(outPath, outName)

        except:

            ARCPY.AddIDMessage("ERROR", 541)

            raise SystemExit()

           

        #### Add Result Fields ####

        fn = UTILS.getFieldNames(outputFieldNames, outPath)

        neighborhoodSizeName, heightMeanName = fn

        outputFields = [neighborhoodSizeName, heightMeanName]

        for field in outputFields:

            UTILS.addEmptyField(outputTable, field, "DOUBLE")

        #### Create Insert Cursor ####

        try:

            insert = DA.InsertCursor(outputTable, outputFields)

        except:

            ARCPY.AddIDMessage("ERROR", 204)

            raise SystemExit()

        #### Add Rows to Output Table ####

        for inc in xrange(self.nIncrements):

            distVal = self.cutoffs[inc]

            ldVal = self.ld[inc]

            rowResult = [distVal, ldVal]

            insert.insertRow(rowResult)

        #### Clean Up ####

        del insert

   

        #### Make Table Visable in TOC if *.dbf Had To Be Added ####

        if dbf:

            ARCPY.SetParameterAsText(1, outputTable)

        ## Save Results ##

        if saveIt:

            myRaster = ARCPY.NumPyArrayToRaster(self.myArray, self.lower_left_corner,

                                                self.xSize, self.ySize,

                                                value_to_nodata = NODATA)

            outPath = OS.path.dirname(self.inRaster)

            outName = "relief" + str(int(self.cutoffs[0])) + ".tif"

            outRaster = OS.path.join(outPath, outName)

            if OS.path.exists(outRaster):

                del outRaster

            else:

                myRaster.save(outRaster)

            ARCPY.AddMessage(myRaster)

           

if __name__ == "__main__":                                                                                                                    第三部分

    setupCostPath()                                                                                                                                     整个程序从这里开始运行，即运行setupCostPath()