栅格数据元数据生成

processing/tools/Raster/test.py

@@ -0,0 +1,140 @@
+from PyQt5.QtCore import QVariant
+
+from qgis.core import *
+from qgis.utils import iface
+from qgis.core import QgsExpression, QgsFeatureRequest, QgsVectorLayer
+
+# 获取渔网图层（假设图层已经加载）
+shp_path = "E:\\projects\\遥感技术部需求\\批量生成元数据\\make\\template\\allfishnet.shp"
+layer = QgsVectorLayer(shp_path, "Fish Layer", "ogr")
+
+# 获取网格宽度和高度（假设是规则网格）
+extent = layer.extent()
+
+grid_width = 100  # 替换为实际网格宽度
+grid_height = 100  # 替换为实际网格高度
+
+# 确保图层有至少一个要素
+if layer.featureCount() > 0:
+    # 获取第一个要素
+    feature = next(layer.getFeatures())
+    geom = feature.geometry()
+    # 获取要素的边界矩形（bounding box）
+    bounds = geom.boundingBox()
+    # 计算宽度和高度
+    grid_width = bounds.width()
+    grid_height = bounds.height()
+    print(f"网格宽度：{grid_width}")
+    print(f"网格高度：{grid_height}")
+else:
+    print("图层没有任何要素！")
+
+# 为图层添加邻居字段
+layer.startEditing()
+for field in ["east", "west", "south", "north"]:
+    if field not in [f.name() for f in layer.fields()]:
+        layer.dataProvider().addAttributes([QgsField(field, QVariant.String)])
+layer.updateFields()
+# 创建空间索引
+spatial_index = QgsSpatialIndex(layer)
+# 设置误差值 过滤掉仅仅是边线重叠的数据
+miss = 0.04
+def cacl(neighbors):
+    if len(neighbors) == 0:
+        return "自由"
+    else:
+        return None
+
+# 遍历图层中的所有要素
+for feature in layer.getFeatures():
+    geom = feature.geometry()
+    extent = geom.boundingBox()
+
+    # 获取网格宽度和高度
+    width = grid_width
+    height = grid_height
+
+    # 构造东西南北方向的矩形
+    west_extent = QgsRectangle(extent.xMinimum() - width + miss, extent.yMinimum() + miss, extent.xMinimum() - miss,
+                               extent.yMaximum() - miss)
+    east_extent = QgsRectangle(extent.xMaximum() + miss, extent.yMinimum() + miss, extent.xMaximum() + width - miss,
+                               extent.yMaximum() - miss)
+    south_extent = QgsRectangle(extent.xMinimum() + miss, extent.yMinimum() - height + miss, extent.xMaximum() - miss,
+                                extent.yMinimum() - miss)
+    north_extent = QgsRectangle(extent.xMinimum() + miss, extent.yMaximum() + miss, extent.xMaximum() - miss,
+                                extent.yMaximum() + height - miss)
+
+    # 查询每个方向是否有其他网格
+    west_neighbors = spatial_index.intersects(west_extent)
+    east_neighbors = spatial_index.intersects(east_extent)
+    south_neighbors = spatial_index.intersects(south_extent)
+    north_neighbors = spatial_index.intersects(north_extent)
+
+    neighbor_geometries = {
+        "north": cacl(north_neighbors),
+        "south": cacl(south_neighbors),
+        "west": cacl(west_neighbors),
+        "east": cacl(east_neighbors)
+    }
+    for direction, value in neighbor_geometries.items():
+        feature[direction] = value
+    layer.updateFeature(feature)
+# 保存编辑
+layer.commitChanges()
+
+# 批量处理每个网格
+# for feature in layer.getFeatures():
+#     geom = feature.geometry()
+#     rectangle = geom.boundingBox()
+#     centroid = geom.centroid().asPoint()
+#     # 定义邻居的目标几何
+#     neighbor_geometries = {
+#         "north": QgsRectangle(
+#             rectangle.xMinimum() + 0,
+#             rectangle.yMinimum() + grid_height,
+#             rectangle.xMaximum() + 0,
+#             rectangle.yMaximum() + grid_height
+#         ),
+#         "south": QgsRectangle(centroid.x() - grid_width / 2, centroid.y() - 1.5 * grid_height,
+#                               centroid.x() + grid_width / 2, centroid.y() - grid_height / 2),
+#         "west": QgsRectangle(centroid.x() - 1.5 * grid_width, centroid.y() - grid_height / 2,
+#                              centroid.x() - grid_width / 2, centroid.y() + grid_height / 2),
+#         "east": QgsRectangle(centroid.x() + grid_width / 2, centroid.y() - grid_height / 2,
+#                              centroid.x() + 1.5 * grid_width, centroid.y() + grid_height / 2)
+#     }
+#     # neighbor_geometries = {
+#     #     "north": QgsRectangle(centroid.x() - grid_width / 2, centroid.y() + grid_height / 2,
+#     #                           centroid.x() + grid_width / 2, centroid.y() + 1.5 * grid_height),
+#     #     "south": QgsRectangle(centroid.x() - grid_width / 2, centroid.y() - 1.5 * grid_height,
+#     #                           centroid.x() + grid_width / 2, centroid.y() - grid_height / 2),
+#     #     "west": QgsRectangle(centroid.x() - 1.5 * grid_width, centroid.y() - grid_height / 2,
+#     #                          centroid.x() - grid_width / 2, centroid.y() + grid_height / 2),
+#     #     "east": QgsRectangle(centroid.x() + grid_width / 2, centroid.y() - grid_height / 2,
+#     #                          centroid.x() + 1.5 * grid_width, centroid.y() + grid_height / 2)
+#     # }
+#     # neighbor_geometries = {
+#     #     "north": QgsRectangle(centroid.x(), centroid.y() + grid_height, centroid.x() + grid_width,
+#     #                            centroid.y() + 2 * grid_height),
+#     #     "south": QgsRectangle(centroid.x(), centroid.y() - 2 * grid_height, centroid.x() + grid_width,
+#     #                              centroid.y() - grid_height),
+#     #     "west": QgsRectangle(centroid.x() - grid_width, centroid.y() - grid_height, centroid.x(),
+#     #                              centroid.y() + grid_height),
+#     #     "east": QgsRectangle(centroid.x() + grid_width, centroid.y() - grid_height, centroid.x() + 2 * grid_width,
+#     #                               centroid.y() + grid_height),
+#     # }
+#     # 查找邻居
+#     neighbors = {}
+#     for direction, neighbor_geom in neighbor_geometries.items():
+#         candidates = spatial_index.intersects(neighbor_geom)
+#         neighbors[direction] = int(
+#             any(layer.getFeature(fid).geometry().intersects(neighbor_geom) for fid in candidates))
+#     # 更新属性
+#     for direction, value in neighbors.items():
+#         if value == 0:
+#             feature[direction] = "自由"
+#         else:
+#             feature[direction] = None
+#     layer.updateFeature(feature)
+# # 保存编辑
+# layer.commitChanges()
+# print("邻居计算完成！")