利用Arcpy和Geopandas進行多規(guī)比對

在ArcGIS中可以很容易進行土規(guī)和城規(guī)的比對,這里,利用Arcpy和Geopandas來進行多規(guī)比對。

土規(guī)
城規(guī)
# encoding: utf-8
# @Author : hanqi

數(shù)據(jù)預處理

導入相關模塊

## 導入相關模塊
import arcpy
from arcpy import env
import numpy as np
import pandas as pd
import geopandas as gpd
import matplotlib.pyplot as plt

防止亂碼

plt.rcParams['font.family'] = ['sans-serif']
plt.rcParams['font.sans-serif'] = ['SimHei']#替換sans-serif字體為黑體
plt.rcParams['axes.unicode_minus'] = False   # 解決坐標軸負數(shù)的負號顯示問題

設置arcpy工作環(huán)境和路徑

env.workspace = r"F:\ArcGIS\多規(guī)對比"
path = "F:\ArcGIS\多規(guī)對比\圖斑對比.gdb"
tg_path = "F:\ArcGIS\多規(guī)對比\圖斑對比.gdb\土規(guī)用地"
cg_path = "F:\ArcGIS\多規(guī)對比\圖斑對比.gdb\城規(guī)用地"

土規(guī)和城規(guī)用地分類

讀取數(shù)據(jù)

land = gpd.GeoDataFrame.from_file(path,layer='土規(guī)用地')
city = gpd.GeoDataFrame.from_file(path,layer='城規(guī)用地')
land.head()
city.head()
處理后的數(shù)據(jù)

arcpy數(shù)據(jù)處理

## 刪除字段
## 因為我做過一遍,所以先刪除
arcpy.DeleteField_management(tg_path, "tg_yongdi")
arcpy.DeleteField_management(cg_path, "cg_yongdi")

## 添加字段
fieldName1 = "tg_yongdi"
fieldLength = 30
 
# Execute AddField twice for two new fields
arcpy.AddField_management(tg_path, fieldName1, "TEXT", field_length=fieldLength)

fieldName2 = "cg_yongdi"
fieldLength = 30

arcpy.AddField_management(cg_path, fieldName2, "TEXT", field_length=fieldLength)
## 字段計算器
## 城規(guī)字段計算器
## 只有一個水域不是建設用地

# Set local variables
inTable = cg_path
fieldName = "cg_yongdi"
expression = "getClass((!cg_fenlei!))"
codeblock = """
def getClass(cg_fenlei):
    if cg_fenlei in "E11":
        return "城規(guī)_非建設用地"
    
    else:
        return "城規(guī)_建設用地"   """
        
# Execute CalculateField 
arcpy.CalculateField_management(inTable, fieldName, expression, "PYTHON_9.3", codeblock)
## 土規(guī)字段計算器

# Set local variables
inTable = tg_path
fieldName = "tg_yongdi"
expression = "getClass((!Layer!))"
codeblock = """
def getClass(Layer):
    if Layer in ("TG-一般農(nóng)田", "TG-基本農(nóng)田", "TG-林地", "TG-水系", "TG-林地"):
        return "土規(guī)_非建設用地"
        
    if Layer is None:
        return None
    
    else:
        return "土規(guī)_建設用地"   """


# Execute CalculateField 
arcpy.CalculateField_management(inTable, fieldName, expression, "PYTHON_9.3", codeblock)

結果就是上圖,土規(guī)和城規(guī)的用地分為兩類了,這里不做演示。

多規(guī)比對

## 聯(lián)合工具

inFeatures = [cg_path, tg_path]
outFeatures = outFeatures = path + "\\多規(guī)比對"

arcpy.Union_analysis (inFeatures, outFeatures, "ALL")
dg_path = "F:\ArcGIS\多歸對比\圖斑對比.gdb\多規(guī)比對"
dg = gpd.GeoDataFrame.from_file(path, layer="多規(guī)比對")
dg.head()
聯(lián)合結果
# Set local variables
## 多規(guī)字段計算器
inTable = dg_path
fieldName = "cg_yongdi"
expression = "getClass3((!cg_yongdi!))"
codeblock = """
def getClass3(cg_yongdi):
    if cg_yongdi is None:
        return "城規(guī)_非建設用地"
    else:
        return cg_yongdi
"""
 
# Execute CalculateField 
arcpy.CalculateField_management(inTable, fieldName, expression, "PYTHON_9.3", codeblock)
## 刪除字段

arcpy.DeleteField_management(dg_path, ["FID_城規(guī)用地","Layer","cg_fenlei","FID_土規(guī)用地"])
## 添加合并字段

fieldName3 = "dg_hebing"
fieldLength = 30
 
# Execute AddField twice for two new fields
arcpy.AddField_management(dg_path, fieldName3, "TEXT", field_length=fieldLength)
## 城規(guī)字段合并
# Set local variables
inTable = dg_path
fieldName = "dg_hebing"
expression = "!tg_yongdi!+'_'+!cg_yongdi!"
# Execute CalculateField 
arcpy.CalculateField_management(inTable, fieldName, expression, "PYTHON_9.3")

dg = gpd.GeoDataFrame.from_file(path, layer="多規(guī)比對")
dg_path = path+'\\多規(guī)比對'
dg
多歸處理后結果

沖突對比

## 添加索引
rec = 0
inTable = dg_path
fieldName = "Identify"
expression = "autoIncrement()"
codeblock = """
def autoIncrement():
    global rec
    pStart = 1
    pInterval = 1
    if (rec == 0):
        rec = pStart
     
    else:
        rec = rec + pInterval
    return rec """
 
# Execute AddField
arcpy.AddField_management(inTable, fieldName, "LONG")
 
# Execute CalculateField 
arcpy.CalculateField_management(inTable, fieldName, expression, "PYTHON_9.3", codeblock)
dg = gpd.GeoDataFrame.from_file(path, layer="多規(guī)比對")
dg.head()
索引
## 按屬性分割


# Set local variables
in_feature_class = dg_path
target_workspace = path
fields = ['dg_hebing']
arcpy.SplitByAttributes_analysis(in_feature_class, target_workspace, fields)
cgfjs_tgjs = gpd.GeoDataFrame.from_file(path, layer="土規(guī)_非建設用地_城規(guī)_建設用地")
cgfjs_tgjs.plot()
土規(guī)_非建設用地_城規(guī)_建設用地
cgjs_tgfjs = gpd.GeoDataFrame.from_file(path, layer="土規(guī)_建設用地_城規(guī)_非建設用地")
cgjs_tgfjs.plot()
土規(guī)_建設用地_城規(guī)_非建設用地

區(qū)劃調(diào)整

土規(guī)非建設用地_城規(guī)建設用地調(diào)整

## 定義路徑
tgfjs_cgjs_path = path+"\\土規(guī)_非建設用地_城規(guī)_建設用地"
# Set local variables
## 土地面積大于10000依城規(guī)

inTable = tgfjs_cgjs_path
fieldName = "gz_1"
expression = "getClass((!Shape_Area!))"
codeblock = """
def getClass(Shape_Area):
    if Shape_Area>10000:
        return "依城規(guī)用地"
        
    else:
        return "依土規(guī)用地"   """


# Execute CalculateField 
arcpy.CalculateField_management(inTable, fieldName, expression, "PYTHON_9.3", codeblock)

土規(guī)_非建設用地_城規(guī)_建設用地屬性表

土規(guī)_非建設用地_城規(guī)_建設用地

土規(guī)建設用地_城規(guī)非建設用地調(diào)整

## 限制1000以內(nèi)不隱藏
pd.set_option('max_columns',1000)
pd.set_option('max_row',1000)
tgjs_cgfjs_path = path+"\\土規(guī)_建設用地_城規(guī)_非建設用地"
tgjs_cgfjs = gpd.GeoDataFrame.from_file(path, layer="土規(guī)_建設用地_城規(guī)_非建設用地")
tgjs_cgfjs.sort_values('Shape_Area',ascending=False)  ##觀察面積
土規(guī)建設_城規(guī)非建設屬性表
# Set local variables
## 土地面積大于10000以土規(guī)

inTable = tgjs_cgfjs
fieldName = "gz_2"
expression = "getClass5((!Shape_Area!))"
codeblock = """
def getClass5(Shape_Area):
    if Shape_Area>10000:
        return "依土規(guī)用地"

    
    else:
        return "依城規(guī)用地"   """


# Execute CalculateField 
arcpy.CalculateField_management(inTable, fieldName, expression, "PYTHON_9.3",   codeblock)
tgjs_cgfjs = gpd.GeoDataFrame.from_file(path, layer="土規(guī)_建設用地_城規(guī)_非建設用地")
tgjs_cgfjs.plot(column='gz_2',legend=True)
土規(guī)建設城規(guī)非建設

數(shù)據(jù)連接


layerName = dg_path
joinTable = cgjs_tgfjs_path
joinField = "Identify"


# Join the feature layer to a table
arcpy.AddJoin_management(layerName, joinField, joinTable, joinField)


layerName = '多規(guī)比對_Layer1'
joinTable = cgfjs_tgjs_path
joinField = "Identify"


# Join the feature layer to a table
arcpy.AddJoin_management(layerName, joinField, joinTable, joinField)

## 導出數(shù)據(jù)
# 輸入
in_features = ['多規(guī)比對_Layer1']
# 輸出
out_location = path
 
# 執(zhí)行 FeatureClassToGeodatabase
arcpy.FeatureClassToGeodatabase_conversion(in_features, out_location)

dg_3 = gpd.GeoDataFrame.from_file(path, layer="多規(guī)比對_Layer3")
dg_3.head()
數(shù)據(jù)連接
## 重命名

# Set local variables
in_data =  "多規(guī)比對_Layer1"
out_data = "多規(guī)比對_最終"

# Execute Rename
arcpy.Rename_management(in_data, out_data)
# Set local variables

inTable = path+"\\多規(guī)對比_最終"
fieldName = "gz_all"
expression = "getClass5(!土規(guī)_非建設用地_城規(guī)_建設用地_gz_1! , !土規(guī)_建設用地_城規(guī)_非建設用地_gz_2!)"  ## 不能用別名
codeblock = """
def getClass6(gz_1,gz_2):
    if gz_1 is not None:
        return gz_1
        
    if gz_2 is not None:
        return gz_2
         
    else:
        return "無沖突"   """

arcpy.AddField_management(inTable, fieldName, "TEXT", 30)

# Execute CalculateField 
arcpy.CalculateField_management(inTable, fieldName, expression, "PYTHON_9.3", codeblock)
dg_zz = gpd.GeoDataFrame.from_file(path, layer="多規(guī)比對_最終")
dg_zz.plot(column='gz_all',legend=True)
沖突比對

用地判定

就是根據(jù)沖突比對來返回想對呀的地塊用地,這里我用字段計算pyhton做不出來,用vb可以,但是可以用pandas來處理。
vb腳本

if [gz_all] = "依土規(guī)用地" then
    value = [多規(guī)比對_tg_yongdi]

if [gz_all] = "依成規(guī)用地" then
    value = [多規(guī)比對_cg_yongdi]

if [gz_all] = "無沖突" then
    value = [多規(guī)比對_cg_yongdi]
end if
value

pandas

## 最終用地
for i in dg_zz.index:
        if dg_zz.at[i,'gz_all']=="無沖突":
            dg_zz["最終用地"].at[i] = dg_zz["多規(guī)比對_tg_yongdi"].at[i]
        if dg_zz.at[i,'gz_all'] == "依城規(guī)用地":
            dg_zz["最終用地"].at[i] = dg_zz["多規(guī)比對_cg_yongdi"].at[i]
        else:
            dg_zz["最終用地"].at[i] = dg_zz["多規(guī)比對_tg_yongdi"].at[i]
dg_zz.head()      
## 分列
dg_zz["用地"] = dg_zz['最終用地'].str.split('_',expand=True)[1]
dg_zz.head()

#dg_zz["用地"] = s.split('_')[1]

結果,這里我進行過很多次計算了


用地結果
## 寫入excel文件

dg_zz.to_excel("成果.xlsx",index=None,encoding="utf8")

然后連接回shp文件就行了,有些ArcGIS不支持xlsx,保存的時候可以注意點,我的可以

字段計算器分列

## 分割字段

# Set local variables

inTable = path + "\\多規(guī)比對_最終"
fieldName = "yongdi"
expression = "getClass6(!zz_yongdi!)"  ## 不能用別名
codeblock = """
def getClass6(zz_yongdi):
    value = zz_yongdi.split('_')[1]
        
    return value   """

arcpy.AddField_management(inTable, fieldName, "TEXT", 30)
# Execute CalculateField 
arcpy.CalculateField_management(inTable, fieldName, expression, "PYTHON_9.3", codeblock)
dg_zz = gpd.GeoDataFrame.from_file(path, layer="多規(guī)比對_最終")
dg_zz.plot(column='zz_yongdi',legend=True)
最終用地比對

成果出圖

fig, ax = plt.subplots(figsize=(10, 10))
ax = dg_zz.plot(ax=ax,column='yongdi',cmap='Set1',legend=True,legend_kwds={
                                                     'loc': 'lower left',
                                                     'title': '圖例',
                                                     'shadow': True})
plt.suptitle('土地利用建設分布圖', fontsize=24) # 添加最高級別標題
plt.tight_layout(pad=4.5) # 調(diào)整不同標題之間間距
plt.grid(True, alpha=0.3)

fig.savefig('土地利用建設分布圖.png', dpi=300)
fig, ax = plt.subplots(figsize=(10, 10))
ax = dg_zz.plot(ax=ax,column='gz_all',cmap='tab10',legend=True,legend_kwds={
                                                     'loc': 'lower left',
                                                     'title': '圖例',
                                                     'shadow': True})
plt.suptitle('規(guī)劃沖突調(diào)整圖', fontsize=24) # 添加最高級別標題
plt.tight_layout(pad=4.5) # 調(diào)整不同標題之間間距
plt.grid(True, alpha=0.3)

fig.savefig('規(guī)劃沖突調(diào)整圖.png', dpi=300)
土地利用建設圖
規(guī)劃沖突調(diào)整圖
?著作權歸作者所有,轉載或內(nèi)容合作請聯(lián)系作者
【社區(qū)內(nèi)容提示】社區(qū)部分內(nèi)容疑似由AI輔助生成,瀏覽時請結合常識與多方信息審慎甄別。
平臺聲明:文章內(nèi)容(如有圖片或視頻亦包括在內(nèi))由作者上傳并發(fā)布,文章內(nèi)容僅代表作者本人觀點,簡書系信息發(fā)布平臺,僅提供信息存儲服務。

相關閱讀更多精彩內(nèi)容

  • 我國現(xiàn)行規(guī)劃體系存在經(jīng)濟社會發(fā)展規(guī)劃和空間規(guī)劃兩大序列,規(guī)劃類型很多(據(jù)不完全統(tǒng)計,具有法定依據(jù)的各類規(guī)劃至少有8...
    ZzxX_b8f7閱讀 1,585評論 0 5
  • 婺源,古徽州一府六縣之一。地處江西省東北部,贛浙皖三省交匯處。 婺源——中國最美的鄉(xiāng)村。一生癡絕處,無夢到徽州。來...
    林可非非非閱讀 6,997評論 2 2
  • 夜鶯2517閱讀 128,174評論 1 9
  • 版本:ios 1.2.1 亮點: 1.app角標可以實時更新天氣溫度或選擇空氣質(zhì)量,建議處女座就不要選了,不然老想...
    我就是沉沉閱讀 7,469評論 1 6
  • 我是一名過去式的高三狗,很可悲,在這三年里我沒有戀愛,看著同齡的小伙伴們一對兒一對兒的,我的心不好受。怎么說呢,高...
    小娘紙閱讀 3,843評論 4 7

友情鏈接更多精彩內(nèi)容