Arcpy实现河流水质动态分段精细化制图

  在博客的“河流水质动态分段精细化制图”中，我们详细介绍了河流水质动态分段精细化制图的流程，下面的内容会通过ArcGIS的脚本工具自动化实现该流程中，根据已经计算了所有河流交汇点的监测数值进行路径创建和动态分段全过程。

输入数据要求：
1、pointFeatureClass：监测站点要素数据（已经计算所有河流交汇点的值）
2、riverPolyline：线状河流数据（与监测站点拓扑关系正确）

一、工具运行界面

  如下图所示，我们为工具定义了多个参数：

二、工具运行结果

  如下图所示，我们预先设置了分级渲染模板，替换数据源即可实现数据更新：

三、脚本实现过程

• 第四步：用监测点与河流交汇点打断河流

# ******************第四步：用监测点与河流交汇点打断河流*******************
#设置处理的容差距离
toleranceValue = "100 Meters"
# 打断的河流要素
outSplitRiver= tempDataPath + "riverpolylinesplit"
arcpy.SplitLineAtPoint_management(riverPolyline, pointFeatureClass, outSplitRiver, toleranceValue)

• 第五步：创建路径

# ******************第五步：创建路径*******************
# 监测站点与打断河流空间相交
outSpatialJoinRiver = tempDataPath + "riverspatialjoinone2one"
arcpy.SpatialJoin_analysis(outSplitRiver, pointFeatureClass, outSpatialJoinRiver, "JOIN_ONE_TO_ONE", "KEEP_ALL", "", "INTERSECT", toleranceValue, "")

# 相交结果增加字段
lengthField = "Length"
arcpy.AddField_management(outSpatialJoinRiver, lengthField, "SHORT", "", "", "", "", "NULLABLE", "NON_REQUIRED", "")
arcpy.CalculateField_management(outSpatialJoinRiver, lengthField, segmentNum, "PYTHON_9.3", "")

# 创建路径
from_measure_field="Length"
outRoute= tempDataPath + "river_createroute"
route_id_field="TARGET_FID"
arcpy.CreateRoutes_lr(outSpatialJoinRiver, route_id_field, outRoute, "ONE_FIELD", from_measure_field, "", "UPPER_LEFT", "1", "0", "IGNORE", "INDEX")

• 第六步：创建路径事件表

# ******************第六步：创建路径事件表*******************
# 路径和监测站点图层空间相交
outRouteSpatialJoin= tempDataPath + "routespatialjoinone2many"
arcpy.SpatialJoin_analysis(outRoute, pointFeatureClass, outRouteSpatialJoin, "JOIN_ONE_TO_MANY", "KEEP_ALL", "", "INTERSECT", toleranceValue, "")

# 获取路径相交结果表中的起点、终点、相交点坐标、相交点监测值等信息
allLineTargetIdList = []
allLineList = []
for row in sorted(arcpy.da.SearchCursor(outRouteSpatialJoin, ["TARGET_FID","SHAPE@LENGTH","SHAPE@","VAL","X","Y"]),key=itemgetter(0),reverse=True):
    targetId = row[0]
    lineGeometry = row[2]
    #获取线段的json格式
    polylineJson = lineGeometry.JSON
    js=json.loads(polylineJson)
    sketchlist =[]
    #读取path
    if ("paths" in js):
        sketchlist = js["paths"]
    else:
        arcpy.AddMessage(targetId)
    if(len(sketchlist)>0):
        #每根线只有一个path
        itemlinePointList = sketchlist[0]
        #起点
        startPoint = itemlinePointList[0]
        startPointX = startPoint[0]
        startPointY = startPoint[1]
        #终点
        endPoint = itemlinePointList[len(itemlinePointList)-1]
        endPointX = endPoint[0]
        endPointY = endPoint[1]
        #相交监测点
        pointX = row[4]
        pointY = row[5]
        pointValue = row[3]
        allLineList.append({"TARGET_FID":targetId,"VAL":pointValue,"X":pointX,"Y":pointY,"startPointX":startPointX,"startPointY":startPointY,"endPointX":endPointX,"endPointY":endPointY})
        allLineTargetIdList.append(targetId)

#根据相交点的坐标判断分段的起始点
resultTableList = []
for targetId in allLineTargetIdList:
    valueStart = 0
    valueEnd = 0
    for lineItem in allLineList:
        if(lineItem["TARGET_FID"]==targetId):
            startPointX = lineItem["startPointX"]
            startPointY = lineItem["startPointY"]
            endPointX = lineItem["endPointX"]
            endPointY = lineItem["endPointY"]
            pointX = lineItem["X"]
            pointY = lineItem["Y"]
            pointValue = lineItem["VAL"]
            #相交点不为空的情况，坐标判断小数点后三位
            if str(pointValue).strip()!='None' and str(pointValue).strip()!='':
                if("%.3f" % startPointX=="%.3f" % pointX and "%.3f" % startPointY=="%.3f" % pointY):
                    valueStart = pointValue
                elif("%.3f" % endPointX=="%.3f" % pointX and "%.3f" % endPointY=="%.3f" % pointY):
                    valueEnd = pointValue
    #每段河流分成100段
    for i in range(100):
        mBeg = i
        mEnd = i+1
        num = (i+1)
        #计算各点的m值 = m值0值点+(m值最大值点-m值0值点)*0.01*k
        resultValue = valueStart +(valueEnd-valueStart)*0.01*num
        obj = (targetId,num,resultValue,valueStart,valueEnd,mBeg,mEnd)
        resultTableList.append(obj)

#路径事件表定义
tableFields = ['TARGET_FID','num','VAL','val_start','val_end','M_beg','M_end']
tabledatatype = numpy.dtype({'names':tableFields,'formats':[numpy.int32,numpy.int32,numpy.float64,numpy.float64,numpy.float64,numpy.float64,numpy.float64]})
tableArr = numpy.array(resultTableList,dtype=tabledatatype)

#如果路径事件表已经存在，需要先删除
routeEventTable = tempDataPath + "routeEventTable"
if arcpy.Exists(routeEventTable):
    arcpy.Delete_management(routeEventTable)
#生成路径事件表
arcpy.da.NumPyArrayToTable(tableArr, routeEventTable)

• 第七步：创建路径事件图层

# ******************第七步：创建路径事件图层*******************
#生成route layer
route_id_field="TARGET_FID"
out_layer = "riverresulttableEvent"+str(random.randint(1,100))
arcpy.MakeRouteEventLayer_lr(outRoute, route_id_field, routeEventTable, "TARGET_FID LINE M_beg M_end", out_layer, "", "NO_ERROR_FIELD", "NO_ANGLE_FIELD", "NORMAL", "ANGLE", "LEFT", "POINT")
#将route layer保存成要素类
arcpy.CopyFeatures_management(out_layer, resultFeatureClass, "", "0", "0", "0")

  关于河流水质动态分段精细化制图，欢迎大家留言，互相交流学习。

  想了解ArcGIS最新的技术动态和环保最新的应用，请关注微信公众号“环保GIS技术与应用”

转载自：https://blog.csdn.net/liuniu1101/article/details/80667429