50 lines
1.8 KiB
Python
50 lines
1.8 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 = 300 / 3.28084 # converting feet to meters
|
||
|
|
COST_PER_METER = 1.5
|
||
|
|
|
||
|
|
# Locate the existing layers
|
||
|
|
home_points_layer = QgsProject.instance().mapLayersByName('HOME_POINTS')[0]
|
||
|
|
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 home point
|
||
|
|
for home_feature in home_points_layer.getFeatures():
|
||
|
|
home_point = home_feature.geometry().asPoint()
|
||
|
|
|
||
|
|
# Initialize nearest neighbor search
|
||
|
|
nearest_neighbor = None
|
||
|
|
nearest_distance = None
|
||
|
|
|
||
|
|
# Search for the nearest pole
|
||
|
|
for pole_feature in poles_layer.getFeatures():
|
||
|
|
pole_point = pole_feature.geometry().asPoint()
|
||
|
|
distance = great_circle((home_point.y(), home_point.x()), (pole_point.y(), pole_point.x())).meters
|
||
|
|
|
||
|
|
# If the pole is within 300 feet and it's closer than the previous nearest neighbor
|
||
|
|
if distance < SEARCH_RADIUS and (nearest_neighbor is None or distance < nearest_distance):
|
||
|
|
nearest_neighbor = pole_feature
|
||
|
|
nearest_distance = distance
|
||
|
|
|
||
|
|
# If a nearest neighbor was found within 300 feet
|
||
|
|
if nearest_neighbor is not None:
|
||
|
|
# Create a new edge feature
|
||
|
|
edge = QgsFeature()
|
||
|
|
edge.setGeometry(QgsGeometry.fromPolylineXY([home_point, nearest_neighbor.geometry().asPoint()]))
|
||
|
|
|
||
|
|
# Calculate cost
|
||
|
|
cost = nearest_distance * COST_PER_METER
|
||
|
|
|
||
|
|
# Set attributes
|
||
|
|
edge.setAttributes(['Aerial Drop', nearest_distance, cost])
|
||
|
|
pr.addFeature(edge)
|
||
|
|
|
||
|
|
# Update the layer's extent when new features have been added
|
||
|
|
edges_layer.updateExtents()
|
||
|
|
print("Done.")
|