FishNet.py 5.7 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142
  1. # -*- coding: utf-8 -*-
  2. __author__ = 'wanger'
  3. __description__ = '处理渔网数据 判断网格是否为自由 增加标签'
  4. __date__ = '2024-11-25'
  5. __copyright__ = '(C) 2024 by siwei'
  6. __revision__ = '1.0'
  7. import os
  8. from PyQt5.QtCore import QVariant
  9. from qgis.PyQt.QtCore import QCoreApplication
  10. from qgis._core import QgsSpatialIndex, QgsField, QgsRectangle
  11. from qgis.core import (QgsProcessing,
  12. QgsVectorLayer,
  13. QgsFeatureSink,
  14. QgsProcessingException,
  15. QgsProcessingAlgorithm,
  16. QgsProcessingParameterFile,
  17. QgsProcessingParameterFeatureSource,
  18. QgsProcessingParameterFeatureSink)
  19. from qgis import processing
  20. class FishNetProcessingAlgorithm(QgsProcessingAlgorithm):
  21. INPUT = 'INPUT'
  22. def tr(self, string):
  23. return QCoreApplication.translate('Processing', string)
  24. def createInstance(self):
  25. return FishNetProcessingAlgorithm()
  26. def name(self):
  27. return 'makeFishNet'
  28. def displayName(self):
  29. return self.tr('处理网格数据')
  30. def group(self):
  31. return self.tr('栅格数据')
  32. def groupId(self):
  33. return 'metadata'
  34. def shortHelpString(self):
  35. return self.tr("将输入的网格数据进行检查,将接边自由的要素做标记。")
  36. def initAlgorithm(self, config=None):
  37. self.addParameter(QgsProcessingParameterFile(
  38. self.INPUT,
  39. '网格数据',
  40. extension='shp' # 设置文件过滤器,只允许选择shp文件
  41. ))
  42. def cacl(self, neighbors):
  43. if len(neighbors) == 0:
  44. return "自由"
  45. else:
  46. return None
  47. # 执行
  48. def processAlgorithm(self, parameters, context, feedback):
  49. shp_file_path = self.parameterAsString(parameters, self.INPUT, context)
  50. # 检查文件路径是否有效
  51. if not os.path.exists(shp_file_path):
  52. feedback.reportError(f"指定的文件路径 {shp_file_path} 不存在!")
  53. return {}
  54. # 加载Shapefile为QgsVectorLayer
  55. layer = QgsVectorLayer(shp_file_path, "Vector Layer", 'ogr')
  56. layer.startEditing()
  57. grid_width = 100 # 替换为实际网格宽度
  58. grid_height = 100 # 替换为实际网格高度
  59. if layer.featureCount() > 0:
  60. # 获取第一个要素
  61. feature = next(layer.getFeatures())
  62. geom = feature.geometry()
  63. # 获取要素的边界矩形(bounding box)
  64. bounds = geom.boundingBox()
  65. # 计算宽度和高度
  66. grid_width = bounds.width()
  67. grid_height = bounds.height()
  68. print(f"网格宽度:{grid_width}")
  69. print(f"网格高度:{grid_height}")
  70. else:
  71. print("图层没有任何要素!")
  72. return {
  73. "状态": "图层没有任何要素!"
  74. }
  75. # 为图层添加邻居字段
  76. layer.startEditing()
  77. for field in ["east", "west", "south", "north"]:
  78. if field not in [f.name() for f in layer.fields()]:
  79. layer.dataProvider().addAttributes([QgsField(field, QVariant.String)])
  80. layer.updateFields()
  81. # 创建空间索引
  82. spatial_index = QgsSpatialIndex(layer)
  83. # 设置误差值 过滤掉仅仅是边线重叠的数据
  84. tolerance = grid_width if grid_width < grid_height else grid_height
  85. tolerance = tolerance / 5
  86. print(f"网格边界容差:{tolerance}")
  87. # 遍历图层中的所有要素
  88. for feature in layer.getFeatures():
  89. geom = feature.geometry()
  90. extent = geom.boundingBox()
  91. # 获取网格宽度和高度
  92. width = grid_width
  93. height = grid_height
  94. # 构造东西南北方向的矩形
  95. west_extent = QgsRectangle(extent.xMinimum() - width + tolerance, extent.yMinimum() + tolerance,
  96. extent.xMinimum() - tolerance,
  97. extent.yMaximum() - tolerance)
  98. east_extent = QgsRectangle(extent.xMaximum() + tolerance, extent.yMinimum() + tolerance,
  99. extent.xMaximum() + width - tolerance,
  100. extent.yMaximum() - tolerance)
  101. south_extent = QgsRectangle(extent.xMinimum() + tolerance, extent.yMinimum() - height + tolerance,
  102. extent.xMaximum() - tolerance,
  103. extent.yMinimum() - tolerance)
  104. north_extent = QgsRectangle(extent.xMinimum() + tolerance, extent.yMaximum() + tolerance,
  105. extent.xMaximum() - tolerance,
  106. extent.yMaximum() + height - tolerance)
  107. # 查询每个方向是否有其他网格
  108. west_neighbors = spatial_index.intersects(west_extent)
  109. east_neighbors = spatial_index.intersects(east_extent)
  110. south_neighbors = spatial_index.intersects(south_extent)
  111. north_neighbors = spatial_index.intersects(north_extent)
  112. neighbor_geometries = {
  113. "north": self.cacl(north_neighbors),
  114. "south": self.cacl(south_neighbors),
  115. "west": self.cacl(west_neighbors),
  116. "east": self.cacl(east_neighbors)
  117. }
  118. for direction, value in neighbor_geometries.items():
  119. feature[direction] = value
  120. layer.updateFeature(feature)
  121. # 保存编辑
  122. layer.commitChanges()
  123. print("邻居计算完成!")
  124. return {
  125. "状态": "处理成功"
  126. }