55 lines
2.1 KiB
Python
55 lines
2.1 KiB
Python
|
from qgis.PyQt.QtCore import QVariant
|
||
|
|
from qgis.core import (QgsFeature, QgsField, QgsGeometry, QgsPointXY, QgsVectorLayer, QgsProject)
|
||
|
|
from geopy.distance import great_circle
|
||
|
|
|
||
|
|
# Constants
|
||
|
|
SEARCH_RADIUS = 50
|
||
|
|
COST_PER_METER = 19.10
|
||
|
|
|
||
|
|
# Locate the existing layers
|
||
|
|
poles_layer = QgsProject.instance().mapLayersByName('POLES')[0]
|
||
|
|
edges_layer = QgsProject.instance().mapLayersByName('EDGES')[0]
|
||
|
|
|
||
|
|
# Define the data provider
|
||
|
|
pr = edges_layer.dataProvider()
|
||
|
|
|
||
|
|
# Iterate over every pole
|
||
|
|
for pole_feature in poles_layer.getFeatures():
|
||
|
|
pole_point = pole_feature.geometry().asPoint()
|
||
|
|
|
||
|
|
# Initialize nearest neighbor search
|
||
|
|
nearest_vertex = None
|
||
|
|
nearest_distance = None
|
||
|
|
|
||
|
|
# Search for the nearest vertex in underground edges
|
||
|
|
for edge_feature in edges_layer.getFeatures():
|
||
|
|
if edge_feature['type'] == 'Underground':
|
||
|
|
# Extract vertices
|
||
|
|
geom = edge_feature.geometry()
|
||
|
|
vertices = geom.asPolyline() if not geom.isMultipart() else [v for part in geom.asMultiPolyline() for v in part]
|
||
|
|
|
||
|
|
for vertex in vertices:
|
||
|
|
vertex_point = QgsPointXY(vertex[0], vertex[1])
|
||
|
|
distance = great_circle((pole_point.y(), pole_point.x()), (vertex_point.y(), vertex_point.x())).meters
|
||
|
|
|
||
|
|
# If the vertex is within 50 feet and it's closer than the previous nearest neighbor
|
||
|
|
if distance < SEARCH_RADIUS and (nearest_vertex is None or distance < nearest_distance):
|
||
|
|
nearest_vertex = vertex_point
|
||
|
|
nearest_distance = distance
|
||
|
|
|
||
|
|
# If a nearest vertex was found within 50 feet
|
||
|
|
if nearest_vertex is not None:
|
||
|
|
# Create a new edge feature
|
||
|
|
edge = QgsFeature()
|
||
|
|
edge.setGeometry(QgsGeometry.fromPolylineXY([pole_point, nearest_vertex]))
|
||
|
|
|
||
|
|
# Calculate cost
|
||
|
|
cost = nearest_distance * COST_PER_METER
|
||
|
|
|
||
|
|
# Set attributes
|
||
|
|
edge.setAttributes(['Transition', nearest_distance, cost])
|
||
|
|
pr.addFeature(edge)
|
||
|
|
|
||
|
|
# Update the layer's extent when new features have been added
|
||
|
|
edges_layer.updateExtents()
|
||
|
|
print("Done.")
|