profile_interpreter.py (23788B)
1 """Profile Interpreter — proof of concept. 2 3 Captures mouse clicks on the QGIS Elevation Profile canvas and places PointZ 4 features at the corresponding along-section position and elevation. This lets 5 you digitize interpretation points (geological boundaries, borehole picks, 6 geophysical layer tops, ...) directly on a profile view. 7 8 The conversion path is: 9 10 click position (canvas pixels, snapped when profile snapping is on) 11 -> QgsElevationProfileCanvas.canvasPointToPlotPoint() 12 -> QgsProfilePoint(distance_along_section, elevation) 13 -> interpolate XY along the profile curve at `distance` 14 -> QgsPoint(x, y, z = elevation) -> PointZ feature 15 """ 16 17 import math 18 import os 19 20 from qgis.PyQt.QtCore import QEvent, QObject, QPointF, Qt, QTimer 21 from qgis.PyQt.QtGui import QIcon 22 from qgis.PyQt.QtWidgets import QAction, QToolBar 23 24 try: 25 _LEFT_BUTTON = Qt.MouseButton.LeftButton # PyQt6 / QGIS 4 26 _RIGHT_BUTTON = Qt.MouseButton.RightButton 27 _CHILD_ADDED = QEvent.Type.ChildAdded 28 except AttributeError: 29 _LEFT_BUTTON = Qt.LeftButton # PyQt5 / QGIS 3 30 _RIGHT_BUTTON = Qt.RightButton 31 _CHILD_ADDED = QEvent.ChildAdded 32 33 from qgis.core import ( 34 Qgis, 35 QgsApplication, 36 QgsCoordinateTransform, 37 QgsFeature, 38 QgsField, 39 QgsGeometry, 40 QgsMessageLog, 41 QgsPoint, 42 QgsPointXY, 43 QgsProfilePoint, 44 QgsProject, 45 QgsVectorDataProvider, 46 QgsVectorLayer, 47 QgsWkbTypes, 48 ) 49 from qgis.PyQt.QtCore import QMetaType 50 from qgis.gui import QgsElevationProfileCanvas, QgsPlotTool 51 52 try: 53 from qgis.gui import QgsPlotToolPan 54 except ImportError: # older builds expose it under a different path 55 QgsPlotToolPan = None 56 57 try: 58 _POINT_GEOMETRY = QgsWkbTypes.GeometryType.PointGeometry # PyQt6 / QGIS 4 59 _ADD_FEATURES_CAP = QgsVectorDataProvider.Capability.AddFeatures 60 _DELETE_FEATURES_CAP = QgsVectorDataProvider.Capability.DeleteFeatures 61 _CHANGE_GEOMETRIES_CAP = QgsVectorDataProvider.Capability.ChangeGeometries 62 _CHANGE_ATTRIBUTES_CAP = QgsVectorDataProvider.Capability.ChangeAttributeValues 63 except AttributeError: 64 _POINT_GEOMETRY = QgsWkbTypes.PointGeometry # PyQt5 / QGIS 3 65 _ADD_FEATURES_CAP = QgsVectorDataProvider.AddFeatures 66 _DELETE_FEATURES_CAP = QgsVectorDataProvider.DeleteFeatures 67 _CHANGE_GEOMETRIES_CAP = QgsVectorDataProvider.ChangeGeometries 68 _CHANGE_ATTRIBUTES_CAP = QgsVectorDataProvider.ChangeAttributeValues 69 70 71 MENU = 'Profile Interpreter' 72 LAYER_NAME = 'Profile interpretations' 73 HIT_TOLERANCE_PX = 10.0 74 75 76 def _is_suitable_target(layer): 77 """True if layer is a writable PointZ vector layer.""" 78 if layer is None or not isinstance(layer, QgsVectorLayer): 79 return False 80 if layer.geometryType() != _POINT_GEOMETRY: 81 return False 82 if not QgsWkbTypes.hasZ(layer.wkbType()): 83 return False 84 caps = layer.dataProvider().capabilities() 85 return bool(caps & _ADD_FEATURES_CAP) 86 87 88 def _point_attributes(field_names, feature_id, distance, elevation): 89 """Return {field: value} for only the names present in field_names.""" 90 values = {'id': feature_id, 'distance': distance, 'elevation': elevation} 91 return {k: v for k, v in values.items() if k in field_names} 92 93 94 def _project_xy(xy, src_crs, dst_crs): 95 """Map x/y from src_crs to dst_crs; identity if either CRS is invalid or they are equal.""" 96 if (not src_crs.isValid() or not dst_crs.isValid() or 97 src_crs == dst_crs): 98 return xy.x(), xy.y() 99 ct = QgsCoordinateTransform(src_crs, dst_crs, QgsProject.instance()) 100 p = ct.transform(QgsPointXY(xy.x(), xy.y())) 101 return p.x(), p.y() 102 103 104 def _next_feature_id(layer): 105 """Next integer id for layer: max existing 'id' + 1, or 1. 106 Returns None if the layer has no 'id' field.""" 107 names = {f.name() for f in layer.fields()} 108 if 'id' not in names: 109 return None 110 max_id = 0 111 for feat in layer.getFeatures(): 112 try: 113 val = feat['id'] 114 except (KeyError, IndexError): 115 continue 116 if val is not None and int(val) > max_id: 117 max_id = int(val) 118 return max_id + 1 119 120 121 def _can_delete(layer): 122 """True if features can be removed from layer's provider.""" 123 if layer is None: 124 return False 125 caps = layer.dataProvider().capabilities() 126 return bool(caps & _DELETE_FEATURES_CAP) 127 128 129 def _can_move(layer): 130 """True if feature geometries can be changed in layer's provider.""" 131 if layer is None: 132 return False 133 caps = layer.dataProvider().capabilities() 134 return bool(caps & _CHANGE_GEOMETRIES_CAP) 135 136 137 def _nearest_hit(candidates, click_x, click_y, tolerance_px): 138 """Feature id of the (fid, px, py) candidate nearest the click within tolerance, or None.""" 139 best_fid = None 140 best_dist = tolerance_px 141 for fid, px, py in candidates: 142 dist = math.hypot(px - click_x, py - click_y) 143 if dist <= best_dist: 144 best_fid = fid 145 best_dist = dist 146 return best_fid 147 148 149 def _added_feature_id(added, fallback): 150 """Provider-assigned id of the feature just written, else the fallback feature's id.""" 151 feature = added[0] if added else fallback 152 return feature.id() 153 154 155 def _event_canvas_point(event): 156 """Canvas position of event, snapped to the plot when profile snapping is on.""" 157 snapped = getattr(event, 'snappedPoint', None) 158 point = snapped() if snapped is not None else event.pos() 159 return QPointF(point.x(), point.y()) 160 161 162 def _find_panel_toolbar(canvas, main_window): 163 """First QToolBar in the panel holding canvas, or None. 164 165 Walks up from the canvas and stops at the main window so the plugin never 166 lands its action on one of the main QGIS toolbars. 167 """ 168 widget = canvas.parentWidget() 169 while widget is not None and widget is not main_window: 170 bars = widget.findChildren(QToolBar) 171 if bars: 172 return bars[0] 173 widget = widget.parentWidget() 174 return None 175 176 177 class _ProfilePickTool(QgsPlotTool): 178 """A plot tool that forwards left-button releases to a callback. 179 180 Right-button releases call on_cancel, if given, before default handling. 181 """ 182 183 def __init__(self, canvas, on_pick, on_cancel=None): 184 super().__init__(canvas, 'Interpret on Profile') 185 self._on_pick = on_pick 186 self._on_cancel = on_cancel 187 188 def plotReleaseEvent(self, event): 189 if event.button() == _LEFT_BUTTON: 190 self._on_pick(event) 191 return 192 if event.button() == _RIGHT_BUTTON and self._on_cancel is not None: 193 self._on_cancel() 194 super().plotReleaseEvent(event) 195 196 197 class _PanelWatcher(QObject): 198 """Event filter that calls on_change once per burst of children added to the watched window.""" 199 200 def __init__(self, on_change): 201 super().__init__() 202 self._on_change = on_change 203 self._pending = False 204 205 def eventFilter(self, watched, event): 206 if event.type() == _CHILD_ADDED and not self._pending: 207 self._pending = True 208 QTimer.singleShot(0, self._fire) 209 return False 210 211 def _fire(self): 212 self._pending = False 213 self._on_change() 214 215 216 class _PanelSession: 217 """Interpretation tool, buttons and undo history bound to one Elevation Profile panel.""" 218 219 def __init__(self, plugin, canvas, toolbar): 220 self._plugin = plugin 221 self.iface = plugin.iface 222 self._canvas = canvas 223 self._toolbar = toolbar 224 self._pick_tool = None 225 self._pan_tool = None 226 self._history = [] 227 self._mode = 'add' 228 self._selected = None 229 230 icon_path = os.path.join(os.path.dirname(__file__), 'icon.svg') 231 self._action = QAction(QIcon(icon_path), 'Interpret on Profile', toolbar) 232 self._action.setToolTip('Place interpretation points on this elevation profile') 233 self._action.setCheckable(True) 234 self._action.toggled.connect(self._set_active) 235 236 self._undo_action = QAction( 237 QgsApplication.getThemeIcon('mActionUndo.svg'), 'Undo Point', toolbar) 238 self._undo_action.setToolTip('Remove the last interpretation point') 239 self._undo_action.setEnabled(False) 240 self._undo_action.triggered.connect(self._undo_last_point) 241 242 self._edit_action = QAction( 243 QgsApplication.getThemeIcon('mActionVertexTool.svg'), 'Edit Point', toolbar) 244 self._edit_action.setToolTip( 245 'Click a point to select it, then click its new position') 246 self._edit_action.setCheckable(True) 247 self._edit_action.setEnabled(False) 248 self._edit_action.toggled.connect(lambda on: self._set_mode('edit', on)) 249 250 self._delete_action = QAction( 251 QgsApplication.getThemeIcon('mActionDeleteSelectedFeatures.svg'), 252 'Delete Point', toolbar) 253 self._delete_action.setToolTip('Click a point to delete it') 254 self._delete_action.setCheckable(True) 255 self._delete_action.setEnabled(False) 256 self._delete_action.toggled.connect(lambda on: self._set_mode('delete', on)) 257 258 toolbar.addAction(self._action) 259 260 def remove(self): 261 """Switch the tool off and take this session's buttons off the panel toolbar.""" 262 self._action.setChecked(False) 263 self._toolbar.removeAction(self._action) 264 265 def _set_active(self, active): 266 if not active: 267 self._deactivate() 268 return 269 self._pick_tool = _ProfilePickTool(self._canvas, self._on_pick, self._clear_selection) 270 self._canvas.setTool(self._pick_tool) 271 self._attach_panel_actions() 272 self.iface.messageBar().pushMessage( 273 MENU, 274 'Click on the elevation profile to place interpretation points.', 275 level=Qgis.Info, 276 ) 277 278 def _deactivate(self): 279 self._detach_panel_actions() 280 if self._pick_tool is not None: 281 if QgsPlotToolPan is not None: 282 self._pan_tool = QgsPlotToolPan(self._canvas) 283 else: 284 self._pan_tool = QgsPlotTool(self._canvas, 'Inactive') 285 self._canvas.setTool(self._pan_tool) 286 self._pick_tool = None 287 288 def _panel_actions(self): 289 return [self._undo_action, self._edit_action, self._delete_action] 290 291 def _attach_panel_actions(self): 292 """Enable the mode actions and put them next to Interpret on the panel toolbar.""" 293 for action in (self._edit_action, self._delete_action): 294 action.setEnabled(True) 295 for action in self._panel_actions(): 296 self._toolbar.addAction(action) 297 298 def _detach_panel_actions(self): 299 for action in (self._edit_action, self._delete_action): 300 action.setChecked(False) 301 action.setEnabled(False) 302 self._clear_selection() 303 self._mode = 'add' 304 for action in self._panel_actions(): 305 self._toolbar.removeAction(action) 306 307 def _set_mode(self, mode, checked): 308 """Switch between add, edit and delete; edit and delete are mutually exclusive.""" 309 self._clear_selection() 310 if not checked: 311 if self._mode == mode: 312 self._mode = 'add' 313 return 314 self._mode = mode 315 other = self._delete_action if mode == 'edit' else self._edit_action 316 other.setChecked(False) 317 hint = ('Click a point to select it, then click its new position.' 318 if mode == 'edit' else 'Click a point to delete it.') 319 self.iface.messageBar().pushMessage(MENU, hint, level=Qgis.Info, duration=3) 320 321 def _clear_selection(self): 322 if self._selected is not None: 323 layer, _ = self._selected 324 if QgsProject.instance().mapLayer(layer.id()) is not None: 325 layer.removeSelection() 326 self._selected = None 327 328 def _refresh_undo_state(self): 329 self._undo_action.setEnabled(bool(self._history)) 330 331 def _refresh_layer(self, layer): 332 layer.updateExtents() 333 layer.triggerRepaint() 334 if hasattr(self._canvas, 'refresh'): 335 self._canvas.refresh() 336 337 def _undo_last_point(self): 338 """Delete the most recently placed point, skipping layers since removed.""" 339 project = QgsProject.instance() 340 while self._history: 341 layer, feature_id = self._history.pop() 342 if project.mapLayer(layer.id()) is None: 343 continue 344 if not _can_delete(layer): 345 self.iface.messageBar().pushMessage( 346 MENU, 347 f'Layer "{layer.name()}" does not support deleting features.', 348 level=Qgis.Warning, 349 ) 350 break 351 if not layer.dataProvider().deleteFeatures([feature_id]): 352 self.iface.messageBar().pushMessage( 353 MENU, 354 'Failed to remove the last point.', 355 level=Qgis.Warning, 356 ) 357 break 358 self._refresh_layer(layer) 359 QgsMessageLog.logMessage( 360 f'removed point {feature_id} from "{layer.name()}"', 361 MENU, level=Qgis.Info, 362 ) 363 self.iface.messageBar().pushMessage( 364 MENU, 365 'Removed last interpretation point.', 366 level=Qgis.Success, 367 duration=2, 368 ) 369 break 370 self._refresh_undo_state() 371 372 def _click_to_point(self, event): 373 """(distance, elevation, xy in canvas CRS) under the click, or None off the section.""" 374 canvas = self._canvas 375 profile_point = canvas.canvasPointToPlotPoint(_event_canvas_point(event)) 376 distance = profile_point.distance() 377 elevation = profile_point.elevation() 378 if math.isnan(distance) or math.isnan(elevation): 379 return None 380 381 curve = canvas.profileCurve() 382 if curve is None: 383 return None 384 point_geom = QgsGeometry(curve.clone()).interpolate(distance) 385 if point_geom.isEmpty(): 386 return None 387 return distance, elevation, point_geom.asPoint() 388 389 def _find_point_at(self, event, layer): 390 """Feature id in layer drawn within HIT_TOLERANCE_PX of the click, or None.""" 391 canvas = self._canvas 392 curve = canvas.profileCurve() 393 if curve is None: 394 return None 395 curve_geom = QgsGeometry(curve.clone()) 396 tolerance = canvas.tolerance() if hasattr(canvas, 'tolerance') else 0.0 397 candidates = [] 398 for feature in layer.getFeatures(): 399 point = feature.geometry().constGet() 400 if not isinstance(point, QgsPoint): 401 continue 402 x, y = _project_xy(point, layer.crs(), canvas.crs()) 403 xy_geom = QgsGeometry(QgsPoint(x, y)) 404 if tolerance > 0 and curve_geom.distance(xy_geom) > tolerance: 405 continue 406 distance = curve_geom.lineLocatePoint(xy_geom) 407 pixel = canvas.plotPointToCanvasPoint(QgsProfilePoint(distance, point.z())) 408 candidates.append((feature.id(), pixel.x(), pixel.y())) 409 pos = event.pos() 410 return _nearest_hit(candidates, pos.x(), pos.y(), HIT_TOLERANCE_PX) 411 412 def _on_pick(self, event): 413 if self._mode == 'edit': 414 self._edit_click(event) 415 elif self._mode == 'delete': 416 self._delete_click(event) 417 else: 418 self._add_point(event) 419 420 def _add_point(self, event): 421 picked = self._click_to_point(event) 422 if picked is None: 423 return 424 distance, elevation, xy = picked 425 canvas = self._canvas 426 427 layer = self._plugin._target_layer(canvas.crs()) 428 x, y = _project_xy(xy, canvas.crs(), layer.crs()) 429 feature = QgsFeature(layer.fields()) 430 feature.setGeometry(QgsGeometry(QgsPoint(x, y, elevation))) 431 432 feature_id = _next_feature_id(layer) 433 attrs = _point_attributes( 434 {f.name() for f in layer.fields()}, 435 feature_id, distance, elevation, 436 ) 437 for k, v in attrs.items(): 438 feature[k] = v 439 440 ok, added = layer.dataProvider().addFeatures([feature]) 441 if not ok: 442 self.iface.messageBar().pushMessage( 443 MENU, 444 'Failed to add feature to layer.', 445 level=Qgis.Warning, 446 ) 447 return 448 self._history.append((layer, _added_feature_id(added, feature))) 449 self._refresh_undo_state() 450 self._refresh_layer(layer) 451 452 id_label = f' #{feature_id}' if feature_id is not None else '' 453 QgsMessageLog.logMessage( 454 f'placed point{id_label} at x={x:.3f} y={y:.3f} z={elevation:.3f} ' 455 f'(distance={distance:.3f})', 456 MENU, level=Qgis.Info, 457 ) 458 self.iface.messageBar().pushMessage( 459 MENU, 460 f'Point at x={xy.x():.1f}, y={xy.y():.1f}, ' 461 f'elevation={elevation:.2f}', 462 level=Qgis.Success, 463 duration=2, 464 ) 465 466 def _edit_click(self, event): 467 if self._selected is None: 468 self._select_for_edit(event) 469 else: 470 self._move_selected(event) 471 472 def _select_for_edit(self, event): 473 layer = self._plugin._search_layer() 474 fid = self._find_point_at(event, layer) if layer is not None else None 475 if fid is None: 476 self.iface.messageBar().pushMessage( 477 MENU, 'No point near the click.', level=Qgis.Info, duration=2) 478 return 479 if not _can_move(layer): 480 self.iface.messageBar().pushMessage( 481 MENU, 482 f'Layer "{layer.name()}" does not support moving features.', 483 level=Qgis.Warning, 484 ) 485 return 486 layer.selectByIds([fid]) 487 self._selected = (layer, fid) 488 self.iface.messageBar().pushMessage( 489 MENU, 'Point selected; click its new position.', level=Qgis.Info, duration=2) 490 491 def _move_selected(self, event): 492 layer, fid = self._selected 493 if QgsProject.instance().mapLayer(layer.id()) is None: 494 self._selected = None 495 return 496 picked = self._click_to_point(event) 497 if picked is None: 498 return 499 distance, elevation, xy = picked 500 x, y = _project_xy(xy, self._canvas.crs(), layer.crs()) 501 provider = layer.dataProvider() 502 geometry = QgsGeometry(QgsPoint(x, y, elevation)) 503 if not provider.changeGeometryValues({fid: geometry}): 504 self.iface.messageBar().pushMessage( 505 MENU, 'Failed to move the point.', level=Qgis.Warning) 506 return 507 508 fields = layer.fields() 509 attrs = _point_attributes({f.name() for f in fields}, None, distance, elevation) 510 attrs.pop('id', None) 511 if attrs and provider.capabilities() & _CHANGE_ATTRIBUTES_CAP: 512 provider.changeAttributeValues( 513 {fid: {fields.indexOf(k): v for k, v in attrs.items()}}) 514 515 self._clear_selection() 516 self._refresh_layer(layer) 517 QgsMessageLog.logMessage( 518 f'moved point {fid} in "{layer.name()}" to x={x:.3f} y={y:.3f} ' 519 f'z={elevation:.3f} (distance={distance:.3f})', 520 MENU, level=Qgis.Info, 521 ) 522 self.iface.messageBar().pushMessage( 523 MENU, 524 f'Moved point to x={xy.x():.1f}, y={xy.y():.1f}, ' 525 f'elevation={elevation:.2f}', 526 level=Qgis.Success, 527 duration=2, 528 ) 529 530 def _delete_click(self, event): 531 layer = self._plugin._search_layer() 532 fid = self._find_point_at(event, layer) if layer is not None else None 533 if fid is None: 534 self.iface.messageBar().pushMessage( 535 MENU, 'No point near the click.', level=Qgis.Info, duration=2) 536 return 537 if not _can_delete(layer): 538 self.iface.messageBar().pushMessage( 539 MENU, 540 f'Layer "{layer.name()}" does not support deleting features.', 541 level=Qgis.Warning, 542 ) 543 return 544 if not layer.dataProvider().deleteFeatures([fid]): 545 self.iface.messageBar().pushMessage( 546 MENU, 'Failed to delete the point.', level=Qgis.Warning) 547 return 548 self._plugin._forget_point(layer, fid) 549 self._refresh_layer(layer) 550 QgsMessageLog.logMessage( 551 f'deleted point {fid} from "{layer.name()}"', MENU, level=Qgis.Info) 552 self.iface.messageBar().pushMessage( 553 MENU, 'Deleted interpretation point.', level=Qgis.Success, duration=2) 554 555 556 class ProfileInterpreterPlugin: 557 def __init__(self, iface): 558 self.iface = iface 559 self._layer = None 560 self._sessions = [] 561 self._watcher = None 562 563 def initGui(self): 564 self._watcher = _PanelWatcher(self._scan_panels) 565 self.iface.mainWindow().installEventFilter(self._watcher) 566 self._scan_panels() 567 568 def unload(self): 569 if self._watcher is not None: 570 self.iface.mainWindow().removeEventFilter(self._watcher) 571 self._watcher = None 572 for session in self._sessions: 573 session.remove() 574 self._sessions = [] 575 576 def _scan_panels(self): 577 """Give every profile panel that has a toolbar and no session its own session.""" 578 if self._watcher is None: 579 return 580 main_window = self.iface.mainWindow() 581 known = [session._canvas for session in self._sessions] 582 for canvas in main_window.findChildren(QgsElevationProfileCanvas): 583 if canvas in known: 584 continue 585 toolbar = _find_panel_toolbar(canvas, main_window) 586 if toolbar is None: 587 continue 588 session = _PanelSession(self, canvas, toolbar) 589 canvas.destroyed.connect(lambda *_, s=session: self._drop_session(s)) 590 self._sessions.append(session) 591 592 def _drop_session(self, session): 593 """Forget a session whose panel was closed; its Qt objects are already gone.""" 594 session._clear_selection() 595 if session in self._sessions: 596 self._sessions.remove(session) 597 598 def _forget_point(self, layer, fid): 599 """Remove a deleted point from every session's undo history.""" 600 for session in self._sessions: 601 session._history = [ 602 (lyr, f) for lyr, f in session._history 603 if not (lyr.id() == layer.id() and f == fid) 604 ] 605 session._refresh_undo_state() 606 607 def _target_layer(self, crs): 608 active = self.iface.activeLayer() 609 if _is_suitable_target(active): 610 return active 611 return self._interpretation_layer(crs) 612 613 def _search_layer(self): 614 """The layer points are added to, without creating the fallback layer.""" 615 active = self.iface.activeLayer() 616 if _is_suitable_target(active): 617 return active 618 if self._layer is not None and QgsProject.instance().mapLayer(self._layer.id()): 619 return self._layer 620 return None 621 622 def _interpretation_layer(self, crs): 623 project = QgsProject.instance() 624 if self._layer is not None and project.mapLayer(self._layer.id()): 625 return self._layer 626 627 uri = f'PointZ?crs={crs.authid()}' if crs.isValid() else 'PointZ' 628 layer = QgsVectorLayer(uri, LAYER_NAME, 'memory') 629 provider = layer.dataProvider() 630 provider.addAttributes([ 631 QgsField('id', QMetaType.Type.Int), 632 QgsField('distance', QMetaType.Type.Double), 633 QgsField('elevation', QMetaType.Type.Double), 634 ]) 635 layer.updateFields() 636 project.addMapLayer(layer) 637 self._layer = layer 638 return layer