ランダムポイント作成


面内に乱点打つのはどうしようかと考えたが当てずっぽうに。
マルチパートだと処理時間長いか。
パートごとに面貼って面積比で打つ方がよいかな。

サンプルそのままならサンプル紹介したほうが速いし
サンプルとほぼ同等のコードを書いても意味ない気がするが、pytでそのまま動くものが少ないので置いとく。
やはり何か重大な不具合でもあってまったく提供されていないのかなぁ。

指定レイヤ(面)に乱点打つコード。
試験データ作成用。

import os.path
import arcpy
import csv
import numpy as np

# Python ツールボックスのツール動作のカスタマイズ
# http://resources.arcgis.com/ja/help/main/10.1/index.html#//00150000002m000000
class Toolbox(object):
    def __init__(self):
        """Define the toolbox (the name of the toolbox is the name of the
        .pyt file)."""
        self.label = "MyToolbox"
        self.alias = ""

        # ツールのリスト設定 ( クラス設定 )
        self.tools = [WktCsvTool]

# 乱点出力
class WktCsvTool(object):
    def __init__(self):
        self.label = "ランダムポイント出力"
        self.description = "ランダムポイント出力"
        self.canRunInBackground = False

    def getParameterInfo(self):
        # Python ツールボックスでのパラメータの定義
        # http://resources.arcgis.com/ja/help/main/10.1/index.html#//001500000028000000
        # getParameterの型一覧
        # http://resources.arcgis.com/ja/help/main/10.1/index.html#//001500000035000000
        # パラメータリファレンス
        # http://resources.arcgis.com/ja/help/main/10.1/index.html#//018z00000063000000

        # レイヤ
        param0 = arcpy.Parameter(
        displayName="入力レイヤ",
        name="in_layer",
        datatype="GPFeatureLayer",
        parameterType="Required",
        direction="Input")
        param0.filter.list = ["Polygon"]

        # 点数フィールド
        param1 = arcpy.Parameter(
        displayName="点数規定フィールド",
        name="in_rndPointNumCol",
        datatype="Field",
        parameterType="Optional",
        direction="Input")

        # 整数型のみ設定可能としてパラメータ0に紐づけ
        param1.filter.list = ['Short', 'Long']
        param1.parameterDependencies = [param0.name]

        # 出力点数
        param2 = arcpy.Parameter(
        displayName="出力点数 [フィールド無指定/値なしの場合]",
        name="in_rndPointNum",
        datatype="GPLong",
        parameterType="Required",
        direction="Input")
        param2.value = 100

        # 出力FeatureClass
        param3 = arcpy.Parameter(
        displayName="出力フィーチャクラス",
        name="outFeatures",
        datatype="DEFeatureClass",
        parameterType="Required",
        direction="Output")

        return [param0, param1, param2, param3]

    def isLicensed(self):
        """Set whether tool is licensed to execute."""
        return True

    def updateParameters(self, parameters):
        """if parameters[1].valueAsText is None:
            parameters[2].enabled = True
        else:
            parameters[2].enabled = False
        """
        return


    def updateMessages(self, parameters):
        """Modify the messages created by internal validation for each tool
        parameter.  This method is called after internal validation."""
        return

    def execute(self, parameters, messages):
        # Python ツールボックスでのメッセージの書き込み
        # http://resources.arcgis.com/ja/help/main/10.1/index.html#//001500000036000000
        try:
            layer = parameters[0].value
            cntCol = parameters[1].valueAsText
            pntNum = parameters[2].value
            outFullName = parameters[3].valueAsText

            messages.addMessage(layer.name)
            if not (cntCol is None):
                messages.addMessage(cntCol)

            messages.addMessage(outFullName)

            geomType = "POINT"
            template = None
            outWs = os.path.dirname(outFullName)
            outFc = os.path.basename(outFullName)

            # messages.addMessage("WS:{0} / Name:{1}".format(outWs,outFc))
            has_m = "DISABLED"
            has_z = "DISABLED"
            spRef = None

            #messages.addMessage(layer.dataSource)
            #if layer.supports("DATASOURCE"):
            dataset = layer.dataSource
            spRef = arcpy.Describe(dataset).spatialReference

            arcpy.CreateFeatureclass_management(outWs, outFc, geomType, template, has_m, has_z, spRef)
            arcpy.AddField_management(outFullName, "orgOid", "LONG", 10, "", "", "元のID", "NON_NULLABLE")

            oidCol = arcpy.Describe(layer).OIDFieldName
            outcols = [oidCol,"SHAPE@"]
            if not (cntCol is None):
                outcols.append(cntCol)

            with arcpy.da.SearchCursor(layer, outcols) as cursor, arcpy.da.InsertCursor(outFullName, ["orgOid","SHAPE@"]) as ins:
                for row in cursor:
                    x = 0
                    loopCnt = pntNum
                    if not (cntCol is None) and not (row[2] is None):
                        loopCnt = row[2]

                    while x < loopCnt:
                        px = row[1].extent.XMin + row[1].extent.width * np.random.random()
                        py = row[1].extent.YMin + row[1].extent.height * np.random.random()
                        pt = arcpy.Point(px,py)
                        
                        #注 面内に含まれるまで繰り返される
                        if row[1].contains(pt):
                            x += 1
                            ins.insertRow([ row[0], pt ])
                        
        except Exception as e:
            messages.AddErrorMessage(e.message)
        return
カテゴリー: 開発 タグ: パーマリンク