さて、本日は頂点の間引き処理を行ってみたいと思います。ArcMapの頂点の間引きというジオプロセシングツールを使えばそれを実現することが可能です。頂点の間引き (Generalize)—ヘルプ | ArcGIS for Desktop
ただ、以下のような問題もあります。
- Standardライセンスからのみ使用可能で、Basicまでしか持っていないという場合は使用不可
今回はBasicライセンスしか持っていないという場合でもそれに対応できるPythonスクリプトを作成してみました。
ソースはGistにも公開しています。https://gist.github.com/sanvarie/51bcdd9e1668081bde75
環境
windows7 64bit
ArcGIS10.2.0
Python2.7.3
処理概要
- フィーチャクラスを作成する(間引いたラインを作成するフィーチャクラス)
- 始点終点は間引き処理対象外とする
- 頂点間の角度と距離を計測して間引き対象/対象外を選別する
- 頂点間の角度の差異が1度以下、かつ、頂点間の長さが1m以上を間引き対象とする
対象データ
こういったデータを作成しました。
各頂点はこのようになっています。
サンプルソース
■generalize.py
# -*- coding: utf-8 -*- import arcpy import math class Generalize(object): def __init__(self,workspace,in_feature,out_feature): arcpy.env.workspace = workspace self.workspace = workspace self.in_feature = in_feature self.out_feature = out_feature self.generalize_list = list() self.check_feature_class() def check_feature_class(self): """out_featureに指定したフィーチャクラスが無い場合、作成します""" exist_out_feature = arcpy.ListFeatureClasses(self.out_feature) if len(exist_out_feature) == 0: try: print u"フィーチャクラスを作成します" arcpy.CreateFeatureclass_management(self.workspace, self.out_feature, "POLYLINE") except Exception: u"フィーチャクラスの作成に失敗しました" def set_generalize_list(self): """角度と距離から削除対象の頂点か判定してリストに格納する self.generalize_list --ラインの始点終点などの情報 [0]:OBJECTID [1]:X [2]:Y [3]:始点X [4]:終点X [5]:削除対象区分 0:削除対象,1:始点で削除非対称,2:終点で削除非対称 [6]:ポイントジオメトリ [7]:角度 [8]:距離 """ for i,pt in enumerate(self.generalize_list): #削除対象のポイントの場合 if pt[5] == 0: #角度計算用に対象の頂点とその前後の頂点の角度を取得する angle_first = self.generalize_list[i-1][7] angle_next = pt[7] angle_last = self.generalize_list[i+1][7] #対象の頂点とその前後のポイントの角度の差が1以下、かつ、対象の頂点間の長さが1以上の場合 if math.fabs(math.fabs(angle_first) - math.fabs(angle_next)) <= 1 \ and math.fabs(math.fabs(angle_next) - math.fabs(angle_last)) <= 1 and pt[8] >= 1: self.generalize_list[i].append(1) else: self.generalize_list[i].append(0) else: self.generalize_list[i].append(0) def calc_angle_length(self): """各頂点の角度や長さなどを計算する関数し、計算結果をリストに格納する関数 self.generalize_list --ラインの始点終点などの情報 [0]:OBJECTID [1]:X [2]:Y [3]:始点X [4]:終点X [5]:削除対象区分 0:削除対象,1:始点で削除非対称,2:終点で削除非対称 [6]:ポイントジオメトリ [7]:角度 [8]:距離 [9]:削除対象区分2 0:削除対象外,1:削除 ←[5]はこれを抽出する前段階のもの """ x2,y2,pg = float(),float(),object() generalize_list2 = list() for i,pt in enumerate(self.generalize_list): #始点、終点でもない場合 if pt[5] == 0: #角度を計算 result = math.atan2(pt[2]-y2,pt[1]-x2) #角度をセット self.generalize_list[i].append(math.degrees(result)) #距離をセット self.generalize_list[i].append(pt[6].distanceTo(pg)) #始点の場合 elif pt[5] == 1: #角度を計算 result = math.atan2(pt[2]-self.generalize_list[i+1][2],pt[1]-self.generalize_list[i+1][1]) #角度をセット if math.degrees(result) < 0: self.generalize_list[i].append(math.degrees(result) + 180) else: self.generalize_list[i].append(math.degrees(result) - 180) #距離は0でOK self.generalize_list[i].append(0) #終点の場合 elif pt[5] == 2: #角度を計算 result = math.atan2(pt[2]-generalize_list2[i-1][2],pt[1]-generalize_list2[i-1][1]) #角度をセット if math.degrees(result) < 0: self.generalize_list[i].append(math.degrees(result)) else: self.generalize_list[i].append(math.degrees(result)) #距離は0でOK self.generalize_list[i].append(0) generalize_list2.append(self.generalize_list[i]) x2,y2,pg = pt[1],pt[2],pt[6] #角度と距離から削除対象の頂点か判定 self.set_generalize_list() def get_point_coord(self): """ラインの頂点ごとの情報をリストに格納する""" lines = [row for row in arcpy.da.SearchCursor(self.in_feature, ["OID@","SHAPE@"])] #頂点ごとにOID、X、Y、始点、終点の情報を格納する for line in lines: for parts in line[1]: for part in parts: self.generalize_list.append([line[0],part.X,part.Y,line[1].firstPoint.X,line[1].lastPoint.X]) def calc_first_last(self): """各頂点が始点なのか終点なのかそれ以外なのか(削除対象)を判定する関数 """ for i,pt in enumerate(self.generalize_list): #対象のポイントのXが始点の場合 if pt[1] == pt[3]: #始点で削除非対象 self.generalize_list[i].append(1) elif pt[1] == pt[4]: #終点で削除非対象 self.generalize_list[i].append(2) else: #削除対象 self.generalize_list[i].append(0) def calc_point_distance(self): """各種関数を呼び出し、self.generalize_listに要素を追加していく""" #ラインの頂点ごとの情報を取得 self.get_point_coord() #各頂点の始点終点それ以外を判定 self.calc_first_last() #頂点のジオメトリ情報を取得 self.create_point_geometry() #頂点の角度と長さの計算 self.calc_angle_length() #削除対象の頂点を間引いてラインを作図する self.generalize_line() def generalize_line(self): """頂点を間引いてラインを作図する関数 self.generalize_list --ラインの始点終点などの情報 [0]:OBJECTID [1]:X [2]:Y [3]:始点X [4]:終点X [5]:削除対象区分 0:削除対象,1:始点で削除非対称,2:終点で削除非対称 [6]:ポイントジオメトリ [7]:角度 [8]:距離 [9]:削除対象区分2 0:削除対象外,1:削除 ←[5]はこれを抽出する前段階のもの """ oid = int() array = arcpy.Array() with arcpy.da.InsertCursor(self.out_feature, ["OID@", "SHAPE@"]) as cur: for i,pt in enumerate(self.generalize_list): #OBJECTIDが次に行ったタイミングでラインを作図する。 if pt[0] != oid and i != 0: cur.insertRow((oid,arcpy.Polyline(array))) array = arcpy.Array() array.add(arcpy.Point(pt[1], pt[2])) else: #削除対象の頂点ではない場合のみ頂点の情報を格納していく if pt[9] != 1: array.add(arcpy.Point(pt[1], pt[2])) oid = pt[0] cur.insertRow((oid,arcpy.Polyline(array))) def create_point_geometry(self): """ポイントジオメトリを作成してリストに格納""" point = arcpy.Point() for i,pt in enumerate(self.generalize_list): point.X = pt[1] point.Y = pt[2] self.generalize_list[i].append(arcpy.PointGeometry(point)) if __name__ == '__main__': #対象のgdb workspace = r"C:\ArcPySample\Generalize.gdb" #間引き元のフィーチャクラス名 in_feature = "Line" #間引き先のフィーチャクラス名 out_feature = "LineGeneralize" #インスタンス化 generalize = Generalize(workspace,in_feature,out_feature) #頂点を間引く generalize.calc_point_distance()
結果を確認すると問題なく、頂点が間引きされていることがわかりました。
課題
- データ量が多い場合は少し時間がかかってしまう
- 属性が存在する場合は、間引いた先のフィーチャクラスにはその属性は引き継がれない(generalize_lineメソッドを修正すればいける)
課題は少しありますが、角度や距離をうまく使用すれば、ArcMapの機能よりも柔軟な使用ができそうですね。簡単ではありますが、本日は以上です。