histogram2d


Histogram2d

# -*- coding: cp932 -*-

import arcpy
import numpy as np

class Toolbox(object):
    def __init__(self):
        self.label = "Toolbox"
        self.alias = ""

        # List of tool classes associated with this toolbox
        self.tools = [histogram2d]

class histogram2d(object):
    def __init__(self):
        self.label = "histogram2d"
        self.description = "histogram2d"
        self.canRunInBackground = False

    def getParameterInfo(self):

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

        wParam = arcpy.Parameter(
        displayName="ピクセル幅",
        name="in_width",
        datatype="GPLong",
        parameterType="Required",
        direction="Input")

        checkHparam = arcpy.Parameter(
        displayName="高さ指定",
        name="in_check_height",
        datatype="GPBoolean",
        parameterType="Optional",
        direction="Input")

        hParam = arcpy.Parameter(
        displayName="ピクセル高さ",
        name="in_height",
        datatype="GPLong",
        parameterType="Optional",
        direction="Input")
        

        # 出力ラスタ
        outRaster = arcpy.Parameter(
        displayName="出力ラスタ",
        name="in_outRaster",
        datatype="DERasterDataset",
        parameterType="Required",
        direction="Output")

        return [ featureLayerParam,wParam, checkHparam,hParam,outRaster ]

    def isLicensed(self):
        return True

    def updateParameters(self, parameters):
        if parameters[2].value == True:
            parameters[3].enabled = True
        else:       
            parameters[3].enabled = False
        return

    def updateMessages(self, parameters):
        if parameters[0].value:
          describe = arcpy.Describe(parameters[0].value)
          if describe.shapeType not in ('Point', 'Multipoint'):
            parameters[0].setWarningMessage('入力がポイントでない場合、重心/中点にて計算され意図しない結果となる可能性があります。')
          else:
            parameters[0].clearMessage()
        return

    def execute(self, parameters, messages):
        feature = parameters[0].value
        pxWidth = parameters[1].value
        isChecked = parameters[2].value
        pxHeight = None
        
        outRaster = parameters[4].valueAsText

        try:
            desc = arcpy.Describe(feature)
            extent = desc.extent
            mapWidth = extent.width
            mapHeight = extent.height
            if isChecked:
                pxHeight = parameters[3].value
            else:
                pxHeight = pxWidth * (mapHeight / mapWidth)
            
            x = arcpy.da.FeatureClassToNumPyArray(feature, ["SHAPE@X"])
            y = arcpy.da.FeatureClassToNumPyArray(feature, ["SHAPE@Y"])
            x = x.astype(np.float64)
            y = y.astype(np.float64)

            x = (((x - int(extent.XMin)) * 100) / int(mapWidth)) * pxWidth 
            y = (((y - int(extent.YMin)) * 100) / int(mapHeight)) * pxHeight

            x = x.astype(np.int32)
            y = y.astype(np.int32)
            x = x.astype(np.int32)
            y = y.astype(np.int32)
            #y = y[-1::-1]

            heatmap, xedges, yedges = np.histogram2d(y, x, bins=( pxHeight, pxWidth))
            #heatmap = np.flipud(heatmap)
            myRaster = arcpy.NumPyArrayToRaster(
                heatmap,
                lower_left_corner = extent.lowerLeft,
                x_cell_size = mapWidth / pxWidth,
                y_cell_size = mapHeight / pxHeight,
                value_to_nodata = 0)
            
            #myRaster.save(outRaster)
            arcpy.Flip_management(myRaster,outRaster)
        except Exception as e:
          messages.AddErrorMessage(e.message)
        return
カテゴリー: 開発 タグ: パーマリンク