Arcpy读取SHP

最新推荐文章于 2024-10-12 08:00:00 发布

Chaoying.

最新推荐文章于 2024-10-12 08:00:00 发布

阅读量1.3k

点赞数

分类专栏： Arcpy与ArcGIS 文章标签： GIS 遥感图像处理 arcpy

本文链接：https://blog.csdn.net/TJLCY/article/details/124039922

版权

Arcpy与ArcGIS 专栏收录该内容

12 篇文章 4 订阅

订阅专栏

该代码实现了一个使用ArcPy库进行几何对象操作的功能，计算两个椭圆的相交部分面积并求得它们的相似性。在grid函数中，遍历中心格网的领域格网，通过计算每个相邻格网椭圆的相似性，得到一个相似性矩阵。最后将结果写入sim.txt文件。整个过程涉及到空间分析和文件路径操作。

摘要由CSDN通过智能技术生成

#coding:utf-8_
import arcpy
import os


def sim(shpFile1, shpFile2):
    ds1 = arcpy.CopyFeatures_management(shpFile1, arcpy.Geometry())
    ds2 = arcpy.CopyFeatures_management(shpFile2, arcpy.Geometry())
    e1 = ds1[0]
    e2 = ds2[0]
    # 求两个椭圆的相交部分
    e_intersection = e1.intersect(e2, 4)
    # 相交部分面积
    area_intersection = e_intersection.area
    # 并集面积
    area_union = e1.area + e2.area - area_intersection
    return area_intersection / area_union


def grid(xno, yno, gridno, girdShp, shpDir):
    '''
    :param xno: 中心格网XNO
    :param yno: 中心格网YNO
    :param gridno: 中心格网编号
    :param girdShp: 格网shp
    :param shpDir: 椭圆SHP文件目录
    :return: 中心格网与领域格网对应椭圆的相似性
    '''
    # 根据XNO YNO编号找到中心格网对应的椭圆文件
    if yno < 100:
        shpName1 = 'tazo_{}0{}_DirectionalDistr.shp'.format(xno, yno)
    else:
        shpName1 = 'tazo_{}{}_DirectionalDistr.shp'.format(xno, yno)
    shpFile1 = os.path.join(shpDir, shpName1)
    if not (os.path.exists(shpFile1)):
        return None
    res = []
    # 遍历中心格网的领域格网计算相似性
    for dx in [-1, 0, 1]:
        for dy in [-1, 0, 1]:
            if dx == 0 and dy == 0:
                continue
            xno_ = xno + dx
            yno_ = yno + dy
            if yno_ < 100:
                shpName2 = 'tazo_{}0{}_DirectionalDistr.shp'.format(xno_, yno_)
            else:
                shpName2 = 'tazo_{}{}_DirectionalDistr.shp'.format(xno_, yno_)
            shpFile2 = os.path.join(shpDir, shpName2)
            if not (os.path.exists(shpFile2)):
                continue
            gridno_ = getgrid(xno_, yno_, gridShp)
            simval = sim(shpFile1, shpFile2)
            res.append([gridno, gridno_, simval])
    return res


def getgrid(xno_, yno_, gridShp):
    '''
    :param xno_: 格网XNO
    :param yno_: 格网YNO
    :param gridShp: 格网shp文件
    :return: 根据XNO YNO查询GRID
    '''
    cursor = arcpy.SearchCursor(gridShp, where_clause="XNO='{}' and YNO='{}'".format(xno_, yno_), fields='GRID')
    for row in cursor:
        return row.getValue('GRID')


if __name__ == "__main__":
    gridShp = r'C:\Users\user\Desktop\5\ZhuHai_2000.shp'
    cursor = arcpy.SearchCursor(gridShp, fields='XNO,YNO,GRID')
    res = []
    for row in cursor:
        xno = int(row.getValue('XNO'))
        yno = int(row.getValue('YNO'))
        gridno = row.getValue('GRID')
        gridres = grid(xno, yno, gridno, gridShp, r'C:\Users\user\Desktop\5\shps')
        if gridres:
            res.extend(gridres)

    # 去除重复计算的格网
    res_ = []
    for gridres in res:
        gridres_ = [gridres[1], gridres[0], gridres[2]]
        if gridres_ in res_:
            continue
        else:
            res_.append(gridres)

    with open('sim.txt', mode='w') as f:
        f.write('GRID_1\tGRID_2\tSIM\n')
        for row in res_:
            f.write("{}\t{}\t{}\n".format(row[0],row[1],row[2]))
    print "ok!"