CreateFeature

本文详细介绍了如何使用CreateFeature和Store方法创建地理数据库特征,并通过InsertCursor实现批量插入。同时,文章阐述了在插入新行或特征时启用加载仅模式以提高性能的方法,以及如何利用缓冲区和刷新技术来优化操作流程。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >

1.Using the CreateFeature and Store methods

public static void CreateFeature(IFeatureClass featureClass, IPolyline polyline)
{
    // Build the feature.
    IFeature feature = featureClass.CreateFeature();
    feature.Shape = polyline;

    // Apply the appropriate subtype to the feature.
    ISubtypes subtypes = (ISubtypes)featureClass;
    IRowSubtypes rowSubtypes = (IRowSubtypes)feature;
    if (subtypes.HasSubtype)
    {
        // In this example, the value of 3 represents the polymer vinyl chloride (PVC) subtype.
        rowSubtypes.SubtypeCode = 3;
    }

    // Initialize any default values the feature has.
    rowSubtypes.InitDefaultValues();

    // Update the value on a string field that indicates who installed the feature.
    int contractorFieldIndex = featureClass.FindField("CONTRACTOR");
    feature.set_Value(contractorFieldIndex, "K Johnston");

    // Commit the new feature to the geodatabase.
    feature.Store();
}

2. Using insert cursors

public static void InsertFeaturesUsingCursor(IFeatureClass featureClass, List <
    IGeometry > geometryList)
{
    using(ComReleaser comReleaser = new ComReleaser())
    {
        // Create a feature buffer.
        IFeatureBuffer featureBuffer = featureClass.CreateFeatureBuffer();
        comReleaser.ManageLifetime(featureBuffer);

        // Create an insert cursor.
        IFeatureCursor insertCursor = featureClass.Insert(true);
        comReleaser.ManageLifetime(insertCursor);

        // All of the features to be created are classified as Primary Highways.
        int typeFieldIndex = featureClass.FindField("TYPE");
        featureBuffer.set_Value(typeFieldIndex, "Primary Highway");
        foreach (IGeometry geometry in geometryList)
        {
            // Set the feature buffer's shape and insert it.
            featureBuffer.Shape = geometry;
            insertCursor.InsertFeature(featureBuffer);
        }

        // Flush the buffer to the geodatabase.
        insertCursor.Flush();
    }
}


3. Using load-only mode to insert rows
The  IFeatureClassLoad interface is an optional interface supported by feature classes in ArcSDE, and feature classes and tables in file geodatabases. It can be used to enable "load-only mode," which can improve the performance of data loading.
Only use load-only mode when inserting new rows or features. Do not use load-only mode when editing existing features.
With ArcSDE, putting a feature class in load-only mode disables updating of the spatial index while inserting features. In a file geodatabase, putting a feature class or table in load-only mode disables the updating of spatial and attribute indexes while inserting rows or features. Taking the feature class or table out of load-only mode rebuilds the indexes.(注意更新索引)
While a feature class or table is in load-only mode, other applications cannot work with the data. Place a feature class or table in load-only mode only after acquiring an exclusive schema lock on the feature class via the  ISchemaLock interface.
public static void LoadOnlyModeInsert(IFeatureClass featureClass, List < IGeometry >
    geometryList)
{
    // Cast the feature class to the IFeatureClassLoad interface.
    IFeatureClassLoad featureClassLoad = (IFeatureClassLoad)featureClass;

    // Acquire an exclusive schema lock for the class.
    ISchemaLock schemaLock = (ISchemaLock)featureClass;
    try
    {
        schemaLock.ChangeSchemaLock(esriSchemaLock.esriExclusiveSchemaLock);

        // Enable load-only mode on the feature class.
        featureClassLoad.LoadOnlyMode = true;
        using(ComReleaser comReleaser = new ComReleaser())
        {
            // Create the feature buffer.
            IFeatureBuffer featureBuffer = featureClass.CreateFeatureBuffer();
            comReleaser.ManageLifetime(featureBuffer);

            // Create an insert cursor.
            IFeatureCursor insertCursor = featureClass.Insert(true);
            comReleaser.ManageLifetime(insertCursor);

            // All of the features to be created are classified as Primary Highways.
            int typeFieldIndex = featureClass.FindField("TYPE");
            featureBuffer.set_Value(typeFieldIndex, "Primary Highway");

            foreach (IGeometry geometry in geometryList)
            {
                // Set the feature buffer's shape and insert it.
                featureBuffer.Shape = geometry;
                insertCursor.InsertFeature(featureBuffer);
            }

            // Flush the buffer to the geodatabase.
            insertCursor.Flush();
        }
    }
    catch (Exception)
    {
        // Handle the failure in a way appropriate to the application.
    }
    finally
    {
        // Disable load-only mode on the feature class.
        featureClassLoad.LoadOnlyMode = false;

        // Demote the exclusive schema lock to a shared lock.
        schemaLock.ChangeSchemaLock(esriSchemaLock.esriSharedSchemaLock);
    }
}

4. Using buffering and flushing
When buffering is used with an insert cursor, a call to InsertRow does not write a new row to the database management system (DBMS), but to a client-side buffer. The rows contained in the buffer are written to the DBMS as a batched operation when ICursor.Flush is called, or when the buffer reaches its maximum size. Because a DBMS write can occur on the Flush call or the InsertRow call (if the buffer has reached its maximum size), both of these calls should have appropriate error handling.
Buffering can only be leveraged during an edit session.


# main_window.py(主窗口逻辑) import os import numpy as np from PySide6.QtWidgets import QMainWindow, QFileDialog, QGraphicsScene, QGraphicsView, QMessageBox, QGraphicsPathItem from PySide6.QtGui import QPainterPath, QPen, QBrush, QAction, QTransform, QImage, QPixmap, QColor from PySide6.QtCore import Qt, QRectF, QPointF from osgeo import ogr, gdal from PySide6.QtWidgets import QInputDialog # 新增输入对话框 # 新增自定义图形项类(用于存储属性) class FeatureItem(QGraphicsPathItem): def __init__(self, path, attributes): super().__init__(path) self.attributes = attributes # 存储属性字典 class MainWindow(QMainWindow): def __init__(self): super().__init__() self.setWindowTitle("GIS软件") self.setGeometry(100, 100, 800, 600) ogr.UseExceptions() self.init_ui() self.scene = QGraphicsScene(self) self.graphicsView.setScene(self.scene) # 新增:存储所有几何边界 self.total_bounds = QRectF() def init_ui(self): self.toolBar = self.addToolBar("工具") self.actionOpen_Vector_Data = QAction("打开矢量数据", self) self.toolBar.addAction(self.actionOpen_Vector_Data) # 新增栅格动作 self.actionOpen_Raster_Data = QAction("打开栅格数据", self) self.toolBar.addAction(self.actionOpen_Raster_Data) # 添加到工具栏 # 新增缓冲区分析按钮 self.actionBuffer_Analysis = QAction("缓冲区分析", self) self.toolBar.addAction(self.actionBuffer_Analysis) self.graphicsView = QGraphicsView() self.setCentralWidget(self.graphicsView) # 新增属性查询按钮 self.actionQuery_Attribute = QAction("属性查询", self) self.toolBar.addAction(self.actionQuery_Attribute) self.actionOpen_Vector_Data.triggered.connect(self.open_vector_data) self.actionOpen_Raster_Data.triggered.connect(self.open_raster_data) # 新增连接 self.actionBuffer_Analysis.triggered.connect(self.buffer_analysis) self.actionQuery_Attribute.triggered.connect(self.enable_query_mode) # 新增鼠标点击事件 self.graphicsView.setMouseTracking(True) self.is_query_mode = False # 新增波段组合按钮 self.actionBand_Combination = QAction("波段组合", self) self.toolBar.addAction(self.actionBand_Combination) self.actionBand_Combination.triggered.connect(self.open_band_combination) # 新增栅格裁剪按钮(在init_ui方法末尾添加) self.actionClip_Raster = QAction("栅格裁剪", self) self.toolBar.addAction(self.actionClip_Raster) self.actionClip_Raster.triggered.connect(self.clip_raster) # 新增连接 self.actionBand_Calculation = QAction("波段运算", self) self.toolBar.addAction(self.actionBand_Calculation) self.actionBand_Calculation.triggered.connect(self.band_calculation) # 新增质心绘制按钮(放在init_ui方法中) self.actionDraw_Centroids = QAction("绘制质心", self) self.toolBar.addAction(self.actionDraw_Centroids) self.actionDraw_Centroids.triggered.connect(self.draw_centroids) self.centroid_items = [] # 新增:存储质心图形项 # 新增空间查询按钮(放在init_ui方法中) self.actionSpatial_Query = QAction("空间查询", self) self.toolBar.addAction(self.actionSpatial_Query) self.actionSpatial_Query.triggered.connect(self.enable_spatial_query_mode) self.is_spatial_query_mode = False self.spatial_query_results = [] # 存储查询结果 # 新增保存按钮(放在init_ui方法中) self.actionSave_Vector = QAction("保存矢量数据", self) self.toolBar.addAction(self.actionSave_Vector) self.actionSave_Vector.triggered.connect(self.save_vector_data) # 新增栅格保存按钮(放在init_ui方法中) self.actionSave_Raster = QAction("保存栅格数据", self) self.toolBar.addAction(self.actionSave_Raster) self.actionSave_Raster.triggered.connect(self.save_raster_data) def open_vector_data(self): file_path, _ = QFileDialog.getOpenFileName( self, "打开矢量文件", "", "Shapefile (*.shp);;GeoJSON (*.geojson);;All Files (*)" ) if file_path: self.load_vector_data(file_path) # 新增:自动缩放视图 self.auto_zoom() def load_vector_data(self, file_path): self.scene.clear() self.total_bounds = QRectF() # 重置边界 try: data_source = ogr.Open(file_path, 0) layer = data_source.GetLayer(0) for feature in layer: geom = feature.GetGeometryRef() path = self.geometry_to_qpainterpath(geom) # 更新总边界 if path.boundingRect().isValid(): self.total_bounds = self.total_bounds.united(path.boundingRect()) pen = QPen(Qt.blue, 1) brush = QBrush(Qt.cyan) self.scene.addPath(path, pen, brush) data_source = None except Exception as e: print(f"加载失败: {str(e)}") self.current_vector_path = file_path # 新增这一行 data_source = None def geometry_to_qpainterpath(self, geom): path = QPainterPath() if geom.GetGeometryType() == ogr.wkbPolygon: for ring in range(geom.GetGeometryCount()): linear_ring = geom.GetGeometryRef(ring) points = linear_ring.GetPoints() if points: path.moveTo(points[0][0], points[0][1]) for p in points[1:]: path.lineTo(p[0], p[1]) path.closeSubpath() elif geom.GetGeometryType() == ogr.wkbLineString: points = geom.GetPoints() if points: path.moveTo(points[0][0], points[0][1]) for p in points[1:]: path.lineTo(p[0], p[1]) elif geom.GetGeometryType() == ogr.wkbPoint: x, y = geom.GetX(), geom.GetY() path.addEllipse(x - 2, y - 2, 4, 4) return path def auto_zoom(self): """自动缩放视图到数据范围并放大2倍""" if not self.total_bounds.isValid(): return # 设置场景边界 self.scene.setSceneRect(self.total_bounds) # 获取视图可视区域 view_rect = self.graphicsView.viewport().rect() # 计算缩放比例(自动适应 + 2倍放大) transform = QTransform() transform.scale(2, 2) # 先放大2倍 # 应用缩放并居中 self.graphicsView.setTransform(transform) self.graphicsView.fitInView(self.total_bounds, Qt.KeepAspectRatio) # 新增缓冲区分析方法 def buffer_analysis(self): """执行缓冲区分析""" if not hasattr(self, 'current_vector_path'): QMessageBox.warning(self, "警告", "请先打开矢量数据文件!") return # 获取缓冲距离 distance, ok = QInputDialog.getDouble( self, "缓冲区分析", "输入缓冲距离(单位与数据坐标系一致):", 0.0, 0 ) if not ok: return try: # 重新打开数据源获取几何 data_source = ogr.Open(self.current_vector_path, 0) layer = data_source.GetLayer(0) # 创建缓冲区路径 buffer_path = QPainterPath() pen = QPen(Qt.red, 2, Qt.DashLine) brush = QBrush(QColor(255, 0, 0, 50)) # 半透明红色填充 for feature in layer: geom = feature.GetGeometryRef() buffer_geom = geom.Buffer(distance) path = self.geometry_to_qpainterpath(buffer_geom) buffer_path.addPath(path) # 添加到场景 self.scene.addPath(buffer_path, pen, brush) # 更新视图边界 if buffer_path.boundingRect().isValid(): self.total_bounds = self.total_bounds.united(buffer_path.boundingRect()) self.auto_zoom() data_source = None except Exception as e: QMessageBox.critical(self, "错误", f"缓冲区分析失败: {str(e)}") def load_vector_data(self, file_path): self.scene.clear() self.total_bounds = QRectF() try: data_source = ogr.Open(file_path, 0) layer = data_source.GetLayer(0) # 获取字段定义 layer_defn = layer.GetLayerDefn() field_names = [layer_defn.GetFieldDefn(i).GetName() for i in range(layer_defn.GetFieldCount())] for feature in layer: geom = feature.GetGeometryRef() path = self.geometry_to_qpainterpath(geom) # 创建属性字典 attributes = { "FID": feature.GetFID(), **{name: feature.GetField(name) for name in field_names} } # 使用自定义图形项 item = FeatureItem(path, attributes) item.setPen(QPen(Qt.blue, 1)) item.setBrush(QBrush(Qt.cyan)) self.scene.addItem(item) if path.boundingRect().isValid(): self.total_bounds = self.total_bounds.united(path.boundingRect()) data_source = None except Exception as e: print(f"加载失败: {str(e)}") self.current_vector_path = file_path data_source = None # 新增属性查询方法 def enable_query_mode(self): """启用属性查询模式""" self.is_query_mode = not self.is_query_mode self.actionQuery_Attribute.setText("退出查询" if self.is_query_mode else "属性查询") self.graphicsView.setCursor(Qt.CrossCursor if self.is_query_mode else Qt.ArrowCursor) # 新增鼠标事件处理 def mousePressEvent(self, event): if self.is_query_mode and event.button() == Qt.LeftButton: scene_pos = self.graphicsView.mapToScene(event.pos()) items = self.scene.items(scene_pos, Qt.IntersectsItemShape, Qt.DescendingOrder) for item in items: if isinstance(item, FeatureItem): # 构建属性信息字符串 info = "\n".join([f"{k}: {v}" for k, v in item.attributes.items()]) QMessageBox.information(self, "要素属性", info) return super().mousePressEvent(event) def draw_centroids(self): """独立质心绘制功能""" if not hasattr(self, 'current_vector_path'): QMessageBox.warning(self, "警告", "请先打开矢量数据文件!") return # 清除已有质心 for item in self.centroid_items: self.scene.removeItem(item) self.centroid_items.clear() try: data_source = ogr.Open(self.current_vector_path, 0) layer = data_source.GetLayer(0) for feature in layer: geom = feature.GetGeometryRef() centroid = geom.Centroid() if centroid: # 创建质心图形项 path = QPainterPath() path.addEllipse( QRectF( centroid.GetX() - 0.3, # 修改为0.3像素半径 centroid.GetY() - 0.3, 0.6, 0.6 # 直径0.6像素 ) ) item = self.scene.addPath( path, QPen(Qt.red, 0.1), QBrush(Qt.red) ) self.centroid_items.append(item) data_source = None self.auto_zoom() except Exception as e: QMessageBox.critical(self, "错误", f"质心绘制失败: {str(e)}") # 新增空间查询模式切换方法 def enable_spatial_query_mode(self): """启用空间查询模式""" self.is_spatial_query_mode = not self.is_spatial_query_mode self.actionSpatial_Query.setText("退出空间查询" if self.is_spatial_query_mode else "空间查询") self.graphicsView.setCursor(Qt.CrossCursor if self.is_spatial_query_mode else Qt.ArrowCursor) if not self.is_spatial_query_mode: self.clear_spatial_query_results() # 新增空间查询处理方法 def mousePressEvent(self, event): if self.is_spatial_query_mode and event.button() == Qt.LeftButton: scene_pos = self.graphicsView.mapToScene(event.pos()) items = self.scene.items(scene_pos, Qt.IntersectsItemShape, Qt.DescendingOrder) for item in items: if isinstance(item, FeatureItem): # 获取空间关系选择 relations = ["相交", "包含", "被包含", "接触", "重叠"] relation, ok = QInputDialog.getItem( self, "空间关系选择", "请选择空间关系:", relations, 0, False ) if not ok: return # 执行空间查询 self.perform_spatial_query(item, relation) return super().mousePressEvent(event) # 新增空间查询核心方法 def perform_spatial_query(self, source_item, relation): """执行空间查询并高亮结果""" self.clear_spatial_query_results() try: # 获取源要素几何 source_geom = self.item_to_geometry(source_item) if not source_geom: return # 获取所有要素 all_items = [item for item in self.scene.items() if isinstance(item, FeatureItem)] # 遍历检查空间关系 for target_item in all_items: target_geom = self.item_to_geometry(target_item) if not target_geom: continue # 执行空间关系判断 if relation == "相交" and source_geom.Intersects(target_geom): self.highlight_item(target_item) elif relation == "包含" and source_geom.Contains(target_geom): self.highlight_item(target_item) elif relation == "被包含" and target_geom.Contains(source_geom): self.highlight_item(target_item) elif relation == "接触" and source_geom.Touches(target_geom): self.highlight_item(target_item) elif relation == "重叠" and source_geom.Overlaps(target_geom): self.highlight_item(target_item) except Exception as e: QMessageBox.critical(self, "错误", f"空间查询失败: {str(e)}") # 新增辅助方法 def item_to_geometry(self, item): """将图形项转换为OGR几何对象""" path = item.path() elements = path.toSubpathPolygons(QTransform()) if not elements: return None # 创建多边形几何 geom = ogr.Geometry(ogr.wkbPolygon) ring = ogr.Geometry(ogr.wkbLinearRing) for point in elements[0]: ring.AddPoint(point.x(), point.y()) ring.CloseRings() geom.AddGeometry(ring) return geom def highlight_item(self, item): """高亮显示查询结果""" original_pen = item.pen() highlight_pen = QPen(Qt.yellow, 3) item.setPen(highlight_pen) self.spatial_query_results.append((item, original_pen)) def clear_spatial_query_results(self): """清除查询结果高亮""" for item, original_pen in self.spatial_query_results: item.setPen(original_pen) self.spatial_query_results.clear() def open_raster_data(self): """打开栅格数据文件""" file_path, _ = QFileDialog.getOpenFileName( self, "打开栅格文件", "", "GeoTIFF (*.tif);;JPEG (*.jpg *.jpeg);;PNG (*.png);;All Files (*)" ) if file_path: try: self.load_raster_data(file_path) self.auto_zoom() except Exception as e: QMessageBox.critical(self, "错误", f"加载栅格失败: {str(e)}") def load_raster_data(self, file_path): """加载栅格数据到视图""" # 打开栅格文件(需要用户修改路径的部分) dataset = gdal.Open(file_path) # 相对路径示例:"./data/raster.tif" # 读取第一个波段 band = dataset.GetRasterBand(1) width = dataset.RasterXSize height = dataset.RasterYSize # 转换为numpy数组 data = band.ReadAsArray() # 创建QImage(注意数据类型转换) if data.dtype == np.uint8: format = QImage.Format.Format_Grayscale8 else: format = QImage.Format.Format_ARGB32 q_img = QImage(data.tobytes(), width, height, format) # 创建像素图项 pixmap = QPixmap.fromImage(q_img) raster_item = self.scene.addPixmap(pixmap) # 处理地理坐标(如果存在) geotransform = dataset.GetGeoTransform() if geotransform: # 计算四个角的坐标 x_origin = geotransform[0] y_origin = geotransform[3] pixel_width = geotransform[1] pixel_height = geotransform[5] # 更新场景边界 x_min = x_origin x_max = x_origin + pixel_width * width y_min = y_origin + pixel_height * height y_max = y_origin self.total_bounds = QRectF( QPointF(x_min, y_min), QPointF(x_max, y_max) ) dataset = None # 关闭数据集 def open_band_combination(self): if not hasattr(self, 'current_raster_path'): QMessageBox.warning(self, "警告", "请先打开栅格数据文件!") return # 复用open_raster_data的逻辑 self.open_raster_data() def open_raster_data(self): file_path, _ = QFileDialog.getOpenFileName( self, "打开栅格文件", "", "GeoTIFF (*.tif);;JPEG (*.jpg *.jpeg);;PNG (*.png);;All Files (*)" ) if file_path: try: dataset = gdal.Open(file_path) num_bands = dataset.RasterCount # 获取用户输入的波段组合 red_band, ok1 = QInputDialog.getInt( self, "波段选择", f"红通道波段号 (1-{num_bands}):", 1, 1, num_bands ) green_band, ok2 = QInputDialog.getInt( self, "波段选择", f"绿通道波段号 (1-{num_bands}):", min(2, num_bands), 1, num_bands ) blue_band, ok3 = QInputDialog.getInt( self, "波段选择", f"蓝通道波段号 (1-{num_bands}):", min(3, num_bands), 1, num_bands ) if not (ok1 and ok2 and ok3): return self.load_raster_data(file_path, red_band, green_band, blue_band) self.auto_zoom() self.current_raster_path = file_path # 新增存储当前路径 except Exception as e: QMessageBox.critical(self, "错误", f"加载栅格失败: {str(e)}") def load_raster_data(self, file_path, red_band=1, green_band=2, blue_band=3): """加载栅格数据到视图(支持波段组合)""" dataset = gdal.Open(file_path) width = dataset.RasterXSize height = dataset.RasterYSize # 读取三个波段数据 def read_band(band_num): band = dataset.GetRasterBand(band_num) data = band.ReadAsArray() # 自动拉伸到0-255范围 data_min = data.min() data_max = data.max() return np.clip(((data - data_min) / (data_max - data_min) * 255), 0, 255).astype(np.uint8) # 合并波段 rgb_array = np.dstack([ read_band(red_band), read_band(green_band), read_band(blue_band) ]) # 创建QImage q_img = QImage( rgb_array.data, width, height, 3 * width, # 每像素3字节(RGB) QImage.Format.Format_RGB888 ) # 创建像素图项 pixmap = QPixmap.fromImage(q_img) self.scene.addPixmap(pixmap) # 处理地理坐标(保持原有逻辑) geotransform = dataset.GetGeoTransform() if geotransform: x_origin = geotransform[0] y_origin = geotransform[3] pixel_width = geotransform[1] pixel_height = geotransform[5] x_min = x_origin x_max = x_origin + pixel_width * width y_min = y_origin + pixel_height * height # 计算下边界 y_max = y_origin # 上边界 # 确保坐标顺序正确 if x_min > x_max: x_min, x_max = x_max, x_min if y_min > y_max: y_min, y_max = y_max, y_min self.total_bounds = QRectF(QPointF(x_min, y_min), QPointF(x_max, y_max)) dataset = None # 新增栅格裁剪方法(必须缩进在类内部) def clip_raster(self): """执行栅格裁剪功能""" if not hasattr(self, 'current_raster_path'): QMessageBox.warning(self, "警告", "请先打开栅格数据文件!") return # 选择裁剪矢量文件 vector_path, _ = QFileDialog.getOpenFileName( self, "选择裁剪区域文件", "", "Shapefile (*.shp);;GeoJSON (*.geojson);;All Files (*)" ) if not vector_path: return try: # 获取原始栅格信息 src_ds = gdal.Open(self.current_raster_path) geotransform = src_ds.GetGeoTransform() proj = src_ds.GetProjection() # 获取矢量范围 vector_ds = ogr.Open(vector_path) layer = vector_ds.GetLayer() feature = layer.GetNextFeature() geom = feature.GetGeometryRef() x_min, x_max, y_min, y_max = geom.GetEnvelope() # 创建临时裁剪结果文件 import os # 确保导入os模块 output_path = os.path.splitext(self.current_raster_path)[0] + "_clipped.tif" # 执行裁剪操作 options = gdal.WarpOptions( format='GTiff', outputBounds=[x_min, y_min, x_max, y_max], dstSRS=proj ) gdal.Warp(output_path, src_ds, options=options) # 加载裁剪结果 self.load_raster_data(output_path) self.auto_zoom() # 清理资源 src_ds = None vector_ds = None except Exception as e: QMessageBox.critical(self, "错误", f"栅格裁剪失败: {str(e)}") # 新增波段运算方法 def band_calculation(self): """执行波段运算(示例为NDVI计算)""" if not hasattr(self, 'current_raster_path'): QMessageBox.warning(self, "警告", "请先打开栅格数据文件!") return try: # 获取用户输入参数 red_band, ok1 = QInputDialog.getInt( self, "波段选择", "输入红波段编号 (1-based):", 1, 1, 100 ) nir_band, ok2 = QInputDialog.getInt( self, "波段选择", "输入近红外波段编号 (1-based):", 4, 1, 100 ) if not (ok1 and ok2): return # 读取栅格数据 dataset = gdal.Open(self.current_raster_path) red = dataset.GetRasterBand(red_band).ReadAsArray() nir = dataset.GetRasterBand(nir_band).ReadAsArray() # 执行NDVI计算 ndvi = np.where( (nir + red) == 0, 0, (nir - red) / (nir + red) ).astype(np.float32) # 创建输出文件 output_path, _ = QFileDialog.getSaveFileName( self, "保存结果", "", "GeoTIFF (*.tif)" ) if not output_path: return # 写入结果 driver = gdal.GetDriverByName('GTiff') out_ds = driver.Create( output_path, dataset.RasterXSize, dataset.RasterYSize, 1, gdal.GDT_Float32 ) out_ds.SetGeoTransform(dataset.GetGeoTransform()) out_ds.SetProjection(dataset.GetProjection()) out_ds.GetRasterBand(1).WriteArray(ndvi) out_ds.FlushCache() # 清理资源 dataset = None out_ds = None QMessageBox.information(self, "成功", f"NDVI计算结果已保存至:\n{output_path}") except Exception as e: QMessageBox.critical(self, "错误", f"波段运算失败: {str(e)}") # 新增保存方法 def save_vector_data(self): """保存当前场景中的矢量数据""" if not hasattr(self, 'current_vector_path'): QMessageBox.warning(self, "警告", "没有可保存的矢量数据!") return # 获取保存路径 file_path, _ = QFileDialog.getSaveFileName( self, "保存矢量文件", "", "Shapefile (*.shp);;GeoJSON (*.geojson);;All Files (*)" ) if not file_path: return try: # 获取原始数据源信息 src_ds = ogr.Open(self.current_vector_path) src_layer = src_ds.GetLayer(0) src_defn = src_layer.GetLayerDefn() # 创建目标数据源 driver = ogr.GetDriverByName("ESRI Shapefile" if file_path.endswith(".shp") else "GeoJSON") if os.path.exists(file_path): driver.DeleteDataSource(file_path) dst_ds = driver.CreateDataSource(file_path) # 创建图层(保持与原始数据相同的坐标系) dst_layer = dst_ds.CreateLayer( "saved_features", srs=src_layer.GetSpatialRef(), geom_type=ogr.wkbPolygon ) # 复制字段定义 for i in range(src_defn.GetFieldCount()): field_defn = src_defn.GetFieldDefn(i) dst_layer.CreateField(field_defn) # 遍历场景中的所有要素项 for item in self.scene.items(): if isinstance(item, FeatureItem): # 创建新要素 feature = ogr.Feature(dst_layer.GetLayerDefn()) # 复制属性 for key, value in item.attributes.items(): if key == "FID": continue # FID通常自动生成 if feature.GetFieldIndex(key) != -1: feature.SetField(key, str(value)) # 转换几何 geom = self.item_to_geometry(item) if geom: feature.SetGeometry(geom) dst_layer.CreateFeature(feature) feature = None # 添加缓冲区要素(如果存在) self.save_additional_features(dst_layer, "buffer") # 添加质心要素(如果存在) self.save_additional_features(dst_layer, "centroid") dst_ds = None src_ds = None QMessageBox.information(self, "成功", f"数据已保存至:\n{file_path}") except Exception as e: QMessageBox.critical(self, "错误", f"保存失败: {str(e)}") # 新增辅助保存方法 def save_additional_features(self, layer, feature_type): """保存附加要素(缓冲区/质心)""" items = [] if feature_type == "buffer": items = [item for item in self.scene.items() if item.pen().style() == Qt.DashLine and item.pen().color() == Qt.red] elif feature_type == "centroid": items = self.centroid_items for item in items: geom = self.item_to_geometry(item) if geom: feature = ogr.Feature(layer.GetLayerDefn()) feature.SetGeometry(geom) # 添加类型标识字段 feature.SetField("FEATURE_TYPE", feature_type.upper()) layer.CreateFeature(feature) feature = None # 在item_to_geometry方法中添加点要素支持 def item_to_geometry(self, item): """增强版几何转换(支持点要素)""" path = item.path() elements = path.toSubpathPolygons(QTransform()) if not elements: # 处理点要素 if isinstance(item, QGraphicsPathItem): path = item.path() if path.elementCount() == 1 and path.elementAt(0).isMoveTo(): pt = path.elementAt(0) geom = ogr.Geometry(ogr.wkbPoint) geom.AddPoint(pt.x, pt.y) return geom return None # 原有多边形处理逻辑 geom = ogr.Geometry(ogr.wkbPolygon) ring = ogr.Geometry(ogr.wkbLinearRing) for point in elements[0]: ring.AddPoint(point.x(), point.y()) ring.CloseRings() geom.AddGeometry(ring) return geom # 新增栅格保存方法 def save_raster_data(self): """保存当前显示的栅格数据""" if not hasattr(self, 'current_raster_path'): QMessageBox.warning(self, "警告", "没有可保存的栅格数据!") return # 获取保存路径和格式 file_path, selected_filter = QFileDialog.getSaveFileName( self, "保存栅格文件", "", "GeoTIFF (*.tif);;JPEG (*.jpg *.jpeg);;PNG (*.png);;All Files (*)" ) if not file_path: return try: # 获取当前显示的栅格数据 dataset = gdal.Open(self.current_raster_path) band_count = dataset.RasterCount width = dataset.RasterXSize height = dataset.RasterYSize # 创建输出数据集 driver = gdal.GetDriverByName(self.get_driver_name(file_path)) out_ds = driver.Create( file_path, width, height, band_count, self.get_gdal_datatype(dataset) ) # 复制地理参考信息 out_ds.SetGeoTransform(dataset.GetGeoTransform()) out_ds.SetProjection(dataset.GetProjection()) # 复制所有波段数据 for band_num in range(1, band_count + 1): in_band = dataset.GetRasterBand(band_num) out_band = out_ds.GetRasterBand(band_num) out_band.WriteArray(in_band.ReadAsArray()) out_band.FlushCache() # 清理资源 out_ds = None dataset = None QMessageBox.information(self, "成功", f"栅格数据已保存至:\n{file_path}") except Exception as e: QMessageBox.critical(self, "错误", f"栅格保存失败: {str(e)}") # 新增辅助方法 def get_driver_name(self, file_path): """根据文件扩展名获取GDAL驱动名称""" ext = os.path.splitext(file_path)[1].lower() return { '.tif': 'GTiff', '.jpg': 'JPEG', '.jpeg': 'JPEG', '.png': 'PNG' }.get(ext, 'GTiff') def get_gdal_datatype(self, dataset): """获取与原始数据集匹配的GDAL数据类型""" sample_band = dataset.GetRasterBand(1) return { gdal.GDT_Byte: gdal.GDT_Byte, gdal.GDT_UInt16: gdal.GDT_UInt16, gdal.GDT_Int16: gdal.GDT_Int16, gdal.GDT_UInt32: gdal.GDT_UInt32, gdal.GDT_Int32: gdal.GDT_Int32, gdal.GDT_Float32: gdal.GDT_Float32, gdal.GDT_Float64: gdal.GDT_Float64 }.get(sample_band.DataType, gdal.GDT_Float32) 你看看我的程序哪里出问题了,为什么属性查询我点了之后什么反应都没有啊,也没有弹窗什么的,之前还有弹窗显示我选中的地方的属性呢
07-18
1. 程序功能概述 (先总体介绍程序具有哪些功能,再依次介绍每一项功能。) 1.1 功能xxx概述 1.2 功能xxx概述 2. 功能实现原理及方法 2.1 功能xxx实现原理及方法 有关键的程序设计思路与方法介绍 2.2 功能xxx实现原理及方法 有关键的程序设计思路与方法介绍 3. 结果展示与说明 3.1功能xxx演示结果与说明 有程序运行结果截图与说明 3.2 功能xxx演示结果与说明 有程序运行结果截图与说明 4. 心得体会 上面这个是我的大作业的程序说明文档的模板,请你根据我的程序以及代码写一个程序说明文档,然后我的功能包括打开矢量数据,打开栅格数据,矢量数据的分析功能包括缓冲区分析,属性查询以及绘制质心,空间查询,保存矢量数据,然后栅格数据的功能包括波段组合,栅格裁剪,以及波段运算,还有一个保存栅格数据,然后这个是我的代码# main_window.py(主窗口逻辑) import os import numpy as np from PySide6.QtWidgets import QMainWindow, QFileDialog, QGraphicsScene, QGraphicsView, QMessageBox, QGraphicsPathItem from PySide6.QtGui import QPainterPath, QPen, QBrush, QAction, QTransform, QImage, QPixmap, QColor from PySide6.QtCore import Qt, QRectF, QPointF from osgeo import ogr, gdal from PySide6.QtWidgets import QInputDialog # 新增输入对话框 # 新增自定义图形项类(用于存储属性) class FeatureItem(QGraphicsPathItem): def __init__(self, path, attributes): super().__init__(path) self.attributes = attributes # 存储属性字典 class MainWindow(QMainWindow): def __init__(self): super().__init__() self.setWindowTitle("GIS软件") self.setGeometry(100, 100, 800, 600) ogr.UseExceptions() self.init_ui() self.scene = QGraphicsScene(self) self.graphicsView.setScene(self.scene) # 新增:存储所有几何边界 self.total_bounds = QRectF() def init_ui(self): self.toolBar = self.addToolBar("工具") self.actionOpen_Vector_Data = QAction("打开矢量数据", self) self.toolBar.addAction(self.actionOpen_Vector_Data) # 新增栅格动作 self.actionOpen_Raster_Data = QAction("打开栅格数据", self) self.toolBar.addAction(self.actionOpen_Raster_Data) # 添加到工具栏 # 新增缓冲区分析按钮 self.actionBuffer_Analysis = QAction("缓冲区分析", self) self.toolBar.addAction(self.actionBuffer_Analysis) self.graphicsView = QGraphicsView() self.setCentralWidget(self.graphicsView) # 新增属性查询按钮 self.actionQuery_Attribute = QAction("属性查询", self) self.toolBar.addAction(self.actionQuery_Attribute) self.actionOpen_Vector_Data.triggered.connect(self.open_vector_data) self.actionOpen_Raster_Data.triggered.connect(self.open_raster_data) # 新增连接 self.actionBuffer_Analysis.triggered.connect(self.buffer_analysis) self.actionQuery_Attribute.triggered.connect(self.enable_query_mode) # 新增鼠标点击事件 self.graphicsView.setMouseTracking(True) self.is_query_mode = False # 新增波段组合按钮 self.actionBand_Combination = QAction("波段组合", self) self.toolBar.addAction(self.actionBand_Combination) self.actionBand_Combination.triggered.connect(self.open_band_combination) # 新增栅格裁剪按钮(在init_ui方法末尾添加) self.actionClip_Raster = QAction("栅格裁剪", self) self.toolBar.addAction(self.actionClip_Raster) self.actionClip_Raster.triggered.connect(self.clip_raster) # 新增连接 self.actionBand_Calculation = QAction("波段运算", self) self.toolBar.addAction(self.actionBand_Calculation) self.actionBand_Calculation.triggered.connect(self.band_calculation) # 新增质心绘制按钮(放在init_ui方法中) self.actionDraw_Centroids = QAction("绘制质心", self) self.toolBar.addAction(self.actionDraw_Centroids) self.actionDraw_Centroids.triggered.connect(self.draw_centroids) self.centroid_items = [] # 新增:存储质心图形项 # 新增空间查询按钮(放在init_ui方法中) self.actionSpatial_Query = QAction("空间查询", self) self.toolBar.addAction(self.actionSpatial_Query) self.actionSpatial_Query.triggered.connect(self.enable_spatial_query_mode) self.is_spatial_query_mode = False self.spatial_query_results = [] # 存储查询结果 # 新增保存按钮(放在init_ui方法中) self.actionSave_Vector = QAction("保存矢量数据", self) self.toolBar.addAction(self.actionSave_Vector) self.actionSave_Vector.triggered.connect(self.save_vector_data) # 新增栅格保存按钮(放在init_ui方法中) self.actionSave_Raster = QAction("保存栅格数据", self) self.toolBar.addAction(self.actionSave_Raster) self.actionSave_Raster.triggered.connect(self.save_raster_data) def open_vector_data(self): file_path, _ = QFileDialog.getOpenFileName( self, "打开矢量文件", "", "Shapefile (*.shp);;GeoJSON (*.geojson);;All Files (*)" ) if file_path: self.load_vector_data(file_path) # 新增:自动缩放视图 self.auto_zoom() def load_vector_data(self, file_path): self.scene.clear() self.total_bounds = QRectF() # 重置边界 try: data_source = ogr.Open(file_path, 0) layer = data_source.GetLayer(0) for feature in layer: geom = feature.GetGeometryRef() path = self.geometry_to_qpainterpath(geom) # 更新总边界 if path.boundingRect().isValid(): self.total_bounds = self.total_bounds.united(path.boundingRect()) pen = QPen(Qt.blue, 1) brush = QBrush(Qt.cyan) self.scene.addPath(path, pen, brush) data_source = None except Exception as e: print(f"加载失败: {str(e)}") self.current_vector_path = file_path # 新增这一行 data_source = None def geometry_to_qpainterpath(self, geom): path = QPainterPath() if geom.GetGeometryType() == ogr.wkbPolygon: for ring in range(geom.GetGeometryCount()): linear_ring = geom.GetGeometryRef(ring) points = linear_ring.GetPoints() if points: path.moveTo(points[0][0], points[0][1]) for p in points[1:]: path.lineTo(p[0], p[1]) path.closeSubpath() elif geom.GetGeometryType() == ogr.wkbLineString: points = geom.GetPoints() if points: path.moveTo(points[0][0], points[0][1]) for p in points[1:]: path.lineTo(p[0], p[1]) elif geom.GetGeometryType() == ogr.wkbPoint: x, y = geom.GetX(), geom.GetY() path.addEllipse(x - 2, y - 2, 4, 4) return path def auto_zoom(self): """自动缩放视图到数据范围并放大2倍""" if not self.total_bounds.isValid(): return # 设置场景边界 self.scene.setSceneRect(self.total_bounds) # 获取视图可视区域 view_rect = self.graphicsView.viewport().rect() # 计算缩放比例(自动适应 + 2倍放大) transform = QTransform() transform.scale(2, 2) # 先放大2倍 # 应用缩放并居中 self.graphicsView.setTransform(transform) self.graphicsView.fitInView(self.total_bounds, Qt.KeepAspectRatio) # 新增缓冲区分析方法 def buffer_analysis(self): """执行缓冲区分析""" if not hasattr(self, 'current_vector_path'): QMessageBox.warning(self, "警告", "请先打开矢量数据文件!") return # 获取缓冲距离 distance, ok = QInputDialog.getDouble( self, "缓冲区分析", "输入缓冲距离(单位与数据坐标系一致):", 0.0, 0 ) if not ok: return try: # 重新打开数据源获取几何 data_source = ogr.Open(self.current_vector_path, 0) layer = data_source.GetLayer(0) # 创建缓冲区路径 buffer_path = QPainterPath() pen = QPen(Qt.red, 2, Qt.DashLine) brush = QBrush(QColor(255, 0, 0, 50)) # 半透明红色填充 for feature in layer: geom = feature.GetGeometryRef() buffer_geom = geom.Buffer(distance) path = self.geometry_to_qpainterpath(buffer_geom) buffer_path.addPath(path) # 添加到场景 self.scene.addPath(buffer_path, pen, brush) # 更新视图边界 if buffer_path.boundingRect().isValid(): self.total_bounds = self.total_bounds.united(buffer_path.boundingRect()) self.auto_zoom() data_source = None except Exception as e: QMessageBox.critical(self, "错误", f"缓冲区分析失败: {str(e)}") def load_vector_data(self, file_path): self.scene.clear() self.total_bounds = QRectF() try: data_source = ogr.Open(file_path, 0) layer = data_source.GetLayer(0) # 获取字段定义 layer_defn = layer.GetLayerDefn() field_names = [layer_defn.GetFieldDefn(i).GetName() for i in range(layer_defn.GetFieldCount())] for feature in layer: geom = feature.GetGeometryRef() path = self.geometry_to_qpainterpath(geom) # 创建属性字典 attributes = { "FID": feature.GetFID(), **{name: feature.GetField(name) for name in field_names} } # 使用自定义图形项 item = FeatureItem(path, attributes) item.setPen(QPen(Qt.blue, 1)) item.setBrush(QBrush(Qt.cyan)) self.scene.addItem(item) if path.boundingRect().isValid(): self.total_bounds = self.total_bounds.united(path.boundingRect()) data_source = None except Exception as e: print(f"加载失败: {str(e)}") self.current_vector_path = file_path data_source = None # 新增属性查询方法 def enable_query_mode(self): """启用属性查询模式""" self.is_query_mode = not self.is_query_mode self.actionQuery_Attribute.setText("退出查询" if self.is_query_mode else "属性查询") self.graphicsView.setCursor(Qt.CrossCursor if self.is_query_mode else Qt.ArrowCursor) # 新增鼠标事件处理 def mousePressEvent(self, event): if self.is_query_mode and event.button() == Qt.LeftButton: scene_pos = self.graphicsView.mapToScene(event.pos()) items = self.scene.items(scene_pos, Qt.IntersectsItemShape, Qt.DescendingOrder) for item in items: if isinstance(item, FeatureItem): # 构建属性信息字符串 info = "\n".join([f"{k}: {v}" for k, v in item.attributes.items()]) QMessageBox.information(self, "要素属性", info) return super().mousePressEvent(event) def draw_centroids(self): """独立质心绘制功能""" if not hasattr(self, 'current_vector_path'): QMessageBox.warning(self, "警告", "请先打开矢量数据文件!") return # 清除已有质心 for item in self.centroid_items: self.scene.removeItem(item) self.centroid_items.clear() try: data_source = ogr.Open(self.current_vector_path, 0) layer = data_source.GetLayer(0) for feature in layer: geom = feature.GetGeometryRef() centroid = geom.Centroid() if centroid: # 创建质心图形项 path = QPainterPath() path.addEllipse( QRectF( centroid.GetX() - 0.3, # 修改为0.3像素半径 centroid.GetY() - 0.3, 0.6, 0.6 # 直径0.6像素 ) ) item = self.scene.addPath( path, QPen(Qt.red, 0.1), QBrush(Qt.red) ) self.centroid_items.append(item) data_source = None self.auto_zoom() except Exception as e: QMessageBox.critical(self, "错误", f"质心绘制失败: {str(e)}") # 新增空间查询模式切换方法 def enable_spatial_query_mode(self): """启用空间查询模式""" self.is_spatial_query_mode = not self.is_spatial_query_mode self.actionSpatial_Query.setText("退出空间查询" if self.is_spatial_query_mode else "空间查询") self.graphicsView.setCursor(Qt.CrossCursor if self.is_spatial_query_mode else Qt.ArrowCursor) if not self.is_spatial_query_mode: self.clear_spatial_query_results() # 新增空间查询处理方法 def mousePressEvent(self, event): if self.is_spatial_query_mode and event.button() == Qt.LeftButton: scene_pos = self.graphicsView.mapToScene(event.pos()) items = self.scene.items(scene_pos, Qt.IntersectsItemShape, Qt.DescendingOrder) for item in items: if isinstance(item, FeatureItem): # 获取空间关系选择 relations = ["相交", "包含", "被包含", "接触", "重叠"] relation, ok = QInputDialog.getItem( self, "空间关系选择", "请选择空间关系:", relations, 0, False ) if not ok: return # 执行空间查询 self.perform_spatial_query(item, relation) return super().mousePressEvent(event) # 新增空间查询核心方法 def perform_spatial_query(self, source_item, relation): """执行空间查询并高亮结果""" self.clear_spatial_query_results() try: # 获取源要素几何 source_geom = self.item_to_geometry(source_item) if not source_geom: return # 获取所有要素 all_items = [item for item in self.scene.items() if isinstance(item, FeatureItem)] # 遍历检查空间关系 for target_item in all_items: target_geom = self.item_to_geometry(target_item) if not target_geom: continue # 执行空间关系判断 if relation == "相交" and source_geom.Intersects(target_geom): self.highlight_item(target_item) elif relation == "包含" and source_geom.Contains(target_geom): self.highlight_item(target_item) elif relation == "被包含" and target_geom.Contains(source_geom): self.highlight_item(target_item) elif relation == "接触" and source_geom.Touches(target_geom): self.highlight_item(target_item) elif relation == "重叠" and source_geom.Overlaps(target_geom): self.highlight_item(target_item) except Exception as e: QMessageBox.critical(self, "错误", f"空间查询失败: {str(e)}") # 新增辅助方法 def item_to_geometry(self, item): """将图形项转换为OGR几何对象""" path = item.path() elements = path.toSubpathPolygons(QTransform()) if not elements: return None # 创建多边形几何 geom = ogr.Geometry(ogr.wkbPolygon) ring = ogr.Geometry(ogr.wkbLinearRing) for point in elements[0]: ring.AddPoint(point.x(), point.y()) ring.CloseRings() geom.AddGeometry(ring) return geom def highlight_item(self, item): """高亮显示查询结果""" original_pen = item.pen() highlight_pen = QPen(Qt.yellow, 3) item.setPen(highlight_pen) self.spatial_query_results.append((item, original_pen)) def clear_spatial_query_results(self): """清除查询结果高亮""" for item, original_pen in self.spatial_query_results: item.setPen(original_pen) self.spatial_query_results.clear() def open_raster_data(self): """打开栅格数据文件""" file_path, _ = QFileDialog.getOpenFileName( self, "打开栅格文件", "", "GeoTIFF (*.tif);;JPEG (*.jpg *.jpeg);;PNG (*.png);;All Files (*)" ) if file_path: try: self.load_raster_data(file_path) self.auto_zoom() except Exception as e: QMessageBox.critical(self, "错误", f"加载栅格失败: {str(e)}") def load_raster_data(self, file_path): """加载栅格数据到视图""" # 打开栅格文件(需要用户修改路径的部分) dataset = gdal.Open(file_path) # 相对路径示例:"./data/raster.tif" # 读取第一个波段 band = dataset.GetRasterBand(1) width = dataset.RasterXSize height = dataset.RasterYSize # 转换为numpy数组 data = band.ReadAsArray() # 创建QImage(注意数据类型转换) if data.dtype == np.uint8: format = QImage.Format.Format_Grayscale8 else: format = QImage.Format.Format_ARGB32 q_img = QImage(data.tobytes(), width, height, format) # 创建像素图项 pixmap = QPixmap.fromImage(q_img) raster_item = self.scene.addPixmap(pixmap) # 处理地理坐标(如果存在) geotransform = dataset.GetGeoTransform() if geotransform: # 计算四个角的坐标 x_origin = geotransform[0] y_origin = geotransform[3] pixel_width = geotransform[1] pixel_height = geotransform[5] # 更新场景边界 x_min = x_origin x_max = x_origin + pixel_width * width y_min = y_origin + pixel_height * height y_max = y_origin self.total_bounds = QRectF( QPointF(x_min, y_min), QPointF(x_max, y_max) ) dataset = None # 关闭数据集 def open_band_combination(self): if not hasattr(self, 'current_raster_path'): QMessageBox.warning(self, "警告", "请先打开栅格数据文件!") return # 复用open_raster_data的逻辑 self.open_raster_data() def open_raster_data(self): file_path, _ = QFileDialog.getOpenFileName( self, "打开栅格文件", "", "GeoTIFF (*.tif);;JPEG (*.jpg *.jpeg);;PNG (*.png);;All Files (*)" ) if file_path: try: dataset = gdal.Open(file_path) num_bands = dataset.RasterCount # 获取用户输入的波段组合 red_band, ok1 = QInputDialog.getInt( self, "波段选择", f"红通道波段号 (1-{num_bands}):", 1, 1, num_bands ) green_band, ok2 = QInputDialog.getInt( self, "波段选择", f"绿通道波段号 (1-{num_bands}):", min(2, num_bands), 1, num_bands ) blue_band, ok3 = QInputDialog.getInt( self, "波段选择", f"蓝通道波段号 (1-{num_bands}):", min(3, num_bands), 1, num_bands ) if not (ok1 and ok2 and ok3): return self.load_raster_data(file_path, red_band, green_band, blue_band) self.auto_zoom() self.current_raster_path = file_path # 新增存储当前路径 except Exception as e: QMessageBox.critical(self, "错误", f"加载栅格失败: {str(e)}") def load_raster_data(self, file_path, red_band=1, green_band=2, blue_band=3): """加载栅格数据到视图(支持波段组合)""" dataset = gdal.Open(file_path) width = dataset.RasterXSize height = dataset.RasterYSize # 读取三个波段数据 def read_band(band_num): band = dataset.GetRasterBand(band_num) data = band.ReadAsArray() # 自动拉伸到0-255范围 data_min = data.min() data_max = data.max() return np.clip(((data - data_min) / (data_max - data_min) * 255), 0, 255).astype(np.uint8) # 合并波段 rgb_array = np.dstack([ read_band(red_band), read_band(green_band), read_band(blue_band) ]) # 创建QImage q_img = QImage( rgb_array.data, width, height, 3 * width, # 每像素3字节(RGB) QImage.Format.Format_RGB888 ) # 创建像素图项 pixmap = QPixmap.fromImage(q_img) self.scene.addPixmap(pixmap) # 处理地理坐标(保持原有逻辑) geotransform = dataset.GetGeoTransform() if geotransform: x_origin = geotransform[0] y_origin = geotransform[3] pixel_width = geotransform[1] pixel_height = geotransform[5] x_min = x_origin x_max = x_origin + pixel_width * width y_min = y_origin + pixel_height * height # 计算下边界 y_max = y_origin # 上边界 # 确保坐标顺序正确 if x_min > x_max: x_min, x_max = x_max, x_min if y_min > y_max: y_min, y_max = y_max, y_min self.total_bounds = QRectF(QPointF(x_min, y_min), QPointF(x_max, y_max)) dataset = None # 新增栅格裁剪方法(必须缩进在类内部) def clip_raster(self): """执行栅格裁剪功能""" if not hasattr(self, 'current_raster_path'): QMessageBox.warning(self, "警告", "请先打开栅格数据文件!") return # 选择裁剪矢量文件 vector_path, _ = QFileDialog.getOpenFileName( self, "选择裁剪区域文件", "", "Shapefile (*.shp);;GeoJSON (*.geojson);;All Files (*)" ) if not vector_path: return try: # 获取原始栅格信息 src_ds = gdal.Open(self.current_raster_path) geotransform = src_ds.GetGeoTransform() proj = src_ds.GetProjection() # 获取矢量范围 vector_ds = ogr.Open(vector_path) layer = vector_ds.GetLayer() feature = layer.GetNextFeature() geom = feature.GetGeometryRef() x_min, x_max, y_min, y_max = geom.GetEnvelope() # 创建临时裁剪结果文件 import os # 确保导入os模块 output_path = os.path.splitext(self.current_raster_path)[0] + "_clipped.tif" # 执行裁剪操作 options = gdal.WarpOptions( format='GTiff', outputBounds=[x_min, y_min, x_max, y_max], dstSRS=proj ) gdal.Warp(output_path, src_ds, options=options) # 加载裁剪结果 self.load_raster_data(output_path) self.auto_zoom() # 清理资源 src_ds = None vector_ds = None except Exception as e: QMessageBox.critical(self, "错误", f"栅格裁剪失败: {str(e)}") # 新增波段运算方法 def band_calculation(self): """执行波段运算(示例为NDVI计算)""" if not hasattr(self, 'current_raster_path'): QMessageBox.warning(self, "警告", "请先打开栅格数据文件!") return try: # 获取用户输入参数 red_band, ok1 = QInputDialog.getInt( self, "波段选择", "输入红波段编号 (1-based):", 1, 1, 100 ) nir_band, ok2 = QInputDialog.getInt( self, "波段选择", "输入近红外波段编号 (1-based):", 4, 1, 100 ) if not (ok1 and ok2): return # 读取栅格数据 dataset = gdal.Open(self.current_raster_path) red = dataset.GetRasterBand(red_band).ReadAsArray() nir = dataset.GetRasterBand(nir_band).ReadAsArray() # 执行NDVI计算 ndvi = np.where( (nir + red) == 0, 0, (nir - red) / (nir + red) ).astype(np.float32) # 创建输出文件 output_path, _ = QFileDialog.getSaveFileName( self, "保存结果", "", "GeoTIFF (*.tif)" ) if not output_path: return # 写入结果 driver = gdal.GetDriverByName('GTiff') out_ds = driver.Create( output_path, dataset.RasterXSize, dataset.RasterYSize, 1, gdal.GDT_Float32 ) out_ds.SetGeoTransform(dataset.GetGeoTransform()) out_ds.SetProjection(dataset.GetProjection()) out_ds.GetRasterBand(1).WriteArray(ndvi) out_ds.FlushCache() # 清理资源 dataset = None out_ds = None QMessageBox.information(self, "成功", f"NDVI计算结果已保存至:\n{output_path}") except Exception as e: QMessageBox.critical(self, "错误", f"波段运算失败: {str(e)}") # 新增保存方法 def save_vector_data(self): """保存当前场景中的矢量数据""" if not hasattr(self, 'current_vector_path'): QMessageBox.warning(self, "警告", "没有可保存的矢量数据!") return # 获取保存路径 file_path, _ = QFileDialog.getSaveFileName( self, "保存矢量文件", "", "Shapefile (*.shp);;GeoJSON (*.geojson);;All Files (*)" ) if not file_path: return try: # 获取原始数据源信息 src_ds = ogr.Open(self.current_vector_path) src_layer = src_ds.GetLayer(0) src_defn = src_layer.GetLayerDefn() # 创建目标数据源 driver = ogr.GetDriverByName("ESRI Shapefile" if file_path.endswith(".shp") else "GeoJSON") if os.path.exists(file_path): driver.DeleteDataSource(file_path) dst_ds = driver.CreateDataSource(file_path) # 创建图层(保持与原始数据相同的坐标系) dst_layer = dst_ds.CreateLayer( "saved_features", srs=src_layer.GetSpatialRef(), geom_type=ogr.wkbPolygon ) # 复制字段定义 for i in range(src_defn.GetFieldCount()): field_defn = src_defn.GetFieldDefn(i) dst_layer.CreateField(field_defn) # 遍历场景中的所有要素项 for item in self.scene.items(): if isinstance(item, FeatureItem): # 创建新要素 feature = ogr.Feature(dst_layer.GetLayerDefn()) # 复制属性 for key, value in item.attributes.items(): if key == "FID": continue # FID通常自动生成 if feature.GetFieldIndex(key) != -1: feature.SetField(key, str(value)) # 转换几何 geom = self.item_to_geometry(item) if geom: feature.SetGeometry(geom) dst_layer.CreateFeature(feature) feature = None # 添加缓冲区要素(如果存在) self.save_additional_features(dst_layer, "buffer") # 添加质心要素(如果存在) self.save_additional_features(dst_layer, "centroid") dst_ds = None src_ds = None QMessageBox.information(self, "成功", f"数据已保存至:\n{file_path}") except Exception as e: QMessageBox.critical(self, "错误", f"保存失败: {str(e)}") # 新增辅助保存方法 def save_additional_features(self, layer, feature_type): """保存附加要素(缓冲区/质心)""" items = [] if feature_type == "buffer": items = [item for item in self.scene.items() if item.pen().style() == Qt.DashLine and item.pen().color() == Qt.red] elif feature_type == "centroid": items = self.centroid_items for item in items: geom = self.item_to_geometry(item) if geom: feature = ogr.Feature(layer.GetLayerDefn()) feature.SetGeometry(geom) # 添加类型标识字段 feature.SetField("FEATURE_TYPE", feature_type.upper()) layer.CreateFeature(feature) feature = None # 在item_to_geometry方法中添加点要素支持 def item_to_geometry(self, item): """增强版几何转换(支持点要素)""" path = item.path() elements = path.toSubpathPolygons(QTransform()) if not elements: # 处理点要素 if isinstance(item, QGraphicsPathItem): path = item.path() if path.elementCount() == 1 and path.elementAt(0).isMoveTo(): pt = path.elementAt(0) geom = ogr.Geometry(ogr.wkbPoint) geom.AddPoint(pt.x, pt.y) return geom return None # 原有多边形处理逻辑 geom = ogr.Geometry(ogr.wkbPolygon) ring = ogr.Geometry(ogr.wkbLinearRing) for point in elements[0]: ring.AddPoint(point.x(), point.y()) ring.CloseRings() geom.AddGeometry(ring) return geom # 新增栅格保存方法 def save_raster_data(self): """保存当前显示的栅格数据""" if not hasattr(self, 'current_raster_path'): QMessageBox.warning(self, "警告", "没有可保存的栅格数据!") return # 获取保存路径和格式 file_path, selected_filter = QFileDialog.getSaveFileName( self, "保存栅格文件", "", "GeoTIFF (*.tif);;JPEG (*.jpg *.jpeg);;PNG (*.png);;All Files (*)" ) if not file_path: return try: # 获取当前显示的栅格数据 dataset = gdal.Open(self.current_raster_path) band_count = dataset.RasterCount width = dataset.RasterXSize height = dataset.RasterYSize # 创建输出数据集 driver = gdal.GetDriverByName(self.get_driver_name(file_path)) out_ds = driver.Create( file_path, width, height, band_count, self.get_gdal_datatype(dataset) ) # 复制地理参考信息 out_ds.SetGeoTransform(dataset.GetGeoTransform()) out_ds.SetProjection(dataset.GetProjection()) # 复制所有波段数据 for band_num in range(1, band_count + 1): in_band = dataset.GetRasterBand(band_num) out_band = out_ds.GetRasterBand(band_num) out_band.WriteArray(in_band.ReadAsArray()) out_band.FlushCache() # 清理资源 out_ds = None dataset = None QMessageBox.information(self, "成功", f"栅格数据已保存至:\n{file_path}") except Exception as e: QMessageBox.critical(self, "错误", f"栅格保存失败: {str(e)}") # 新增辅助方法 def get_driver_name(self, file_path): """根据文件扩展名获取GDAL驱动名称""" ext = os.path.splitext(file_path)[1].lower() return { '.tif': 'GTiff', '.jpg': 'JPEG', '.jpeg': 'JPEG', '.png': 'PNG' }.get(ext, 'GTiff') def get_gdal_datatype(self, dataset): """获取与原始数据集匹配的GDAL数据类型""" sample_band = dataset.GetRasterBand(1) return { gdal.GDT_Byte: gdal.GDT_Byte, gdal.GDT_UInt16: gdal.GDT_UInt16, gdal.GDT_Int16: gdal.GDT_Int16, gdal.GDT_UInt32: gdal.GDT_UInt32, gdal.GDT_Int32: gdal.GDT_Int32, gdal.GDT_Float32: gdal.GDT_Float32, gdal.GDT_Float64: gdal.GDT_Float64 }.get(sample_band.DataType, gdal.GDT_Float32)
07-18
这些更是一点儿质心都看不着了,看看代码,问题出在哪# main_window.py(主窗口逻辑) import os import numpy as np from PySide6.QtWidgets import QMainWindow, QFileDialog, QGraphicsScene, QGraphicsView, QMessageBox, QGraphicsPathItem from PySide6.QtGui import QPainterPath, QPen, QBrush, QAction, QTransform, QImage, QPixmap, QColor, QPainter from PySide6.QtCore import Qt, QRectF, QPointF from osgeo import ogr, gdal from PySide6.QtWidgets import QInputDialog # 新增输入对话框 # 新增自定义图形项类(用于存储属性) class FeatureItem(QGraphicsPathItem): def __init__(self, path, attributes): super().__init__(path) self.attributes = attributes # 存储属性字典 class MainWindow(QMainWindow): def __init__(self): super().__init__() self.setWindowTitle("GIS软件") self.setGeometry(100, 100, 800, 600) ogr.UseExceptions() self.init_ui() self.scene = QGraphicsScene(self) self.graphicsView.setScene(self.scene) # 新增:存储所有几何边界 self.total_bounds = QRectF() def init_ui(self): self.toolBar = self.addToolBar("工具") self.actionOpen_Vector_Data = QAction("打开矢量数据", self) self.toolBar.addAction(self.actionOpen_Vector_Data) # 新增栅格动作 self.actionOpen_Raster_Data = QAction("打开栅格数据", self) self.toolBar.addAction(self.actionOpen_Raster_Data) # 添加到工具栏 # 新增缓冲区分析按钮 self.actionBuffer_Analysis = QAction("缓冲区分析", self) self.toolBar.addAction(self.actionBuffer_Analysis) self.graphicsView = QGraphicsView() self.setCentralWidget(self.graphicsView) # 新增属性查询按钮 self.actionQuery_Attribute = QAction("属性查询", self) self.toolBar.addAction(self.actionQuery_Attribute) self.actionOpen_Vector_Data.triggered.connect(self.open_vector_data) self.actionOpen_Raster_Data.triggered.connect(self.open_raster_data) # 新增连接 self.actionBuffer_Analysis.triggered.connect(self.buffer_analysis) self.actionQuery_Attribute.triggered.connect(self.enable_query_mode) # 新增鼠标点击事件 self.graphicsView.setMouseTracking(True) self.is_query_mode = False # 新增波段组合按钮 self.actionBand_Combination = QAction("波段组合", self) self.toolBar.addAction(self.actionBand_Combination) self.actionBand_Combination.triggered.connect(self.open_band_combination) # 新增栅格裁剪按钮(在init_ui方法末尾添加) self.actionClip_Raster = QAction("栅格裁剪", self) self.toolBar.addAction(self.actionClip_Raster) self.actionClip_Raster.triggered.connect(self.clip_raster) # 新增连接 self.actionBand_Calculation = QAction("波段运算", self) self.toolBar.addAction(self.actionBand_Calculation) self.actionBand_Calculation.triggered.connect(self.band_calculation) # 新增质心绘制按钮(放在init_ui方法中) self.actionDraw_Centroids = QAction("绘制质心", self) self.toolBar.addAction(self.actionDraw_Centroids) self.actionDraw_Centroids.triggered.connect(self.draw_centroids) self.centroid_items = [] # 新增:存储质心图形项 # 在init_ui方法中添加 self.graphicsView.setRenderHint(QPainter.Antialiasing) self.graphicsView.setTransformationAnchor(QGraphicsView.AnchorUnderMouse) # 新增空间查询按钮(放在init_ui方法中) self.actionSpatial_Query = QAction("空间查询", self) self.toolBar.addAction(self.actionSpatial_Query) self.actionSpatial_Query.triggered.connect(self.enable_spatial_query_mode) self.is_spatial_query_mode = False self.spatial_query_results = [] # 存储查询结果 # 新增保存按钮(放在init_ui方法中) self.actionSave_Vector = QAction("保存矢量数据", self) self.toolBar.addAction(self.actionSave_Vector) self.actionSave_Vector.triggered.connect(self.save_vector_data) # 新增栅格保存按钮(放在init_ui方法中) self.actionSave_Raster = QAction("保存栅格数据", self) self.toolBar.addAction(self.actionSave_Raster) self.actionSave_Raster.triggered.connect(self.save_raster_data) def open_vector_data(self): file_path, _ = QFileDialog.getOpenFileName( self, "打开矢量文件", "", "Shapefile (*.shp);;GeoJSON (*.geojson);;All Files (*)" ) if file_path: self.load_vector_data(file_path) # 新增:自动缩放视图 self.auto_zoom() def load_vector_data(self, file_path): self.scene.clear() self.total_bounds = QRectF() # 重置边界 try: data_source = ogr.Open(file_path, 0) layer = data_source.GetLayer(0) for feature in layer: geom = feature.GetGeometryRef() path = self.geometry_to_qpainterpath(geom) # 更新总边界 if path.boundingRect().isValid(): self.total_bounds = self.total_bounds.united(path.boundingRect()) pen = QPen(Qt.blue, 1) brush = QBrush(Qt.cyan) self.scene.addPath(path, pen, brush) data_source = None except Exception as e: print(f"加载失败: {str(e)}") self.current_vector_path = file_path # 新增这一行 data_source = None def geometry_to_qpainterpath(self, geom): path = QPainterPath() if geom.GetGeometryType() == ogr.wkbPolygon: for ring in range(geom.GetGeometryCount()): linear_ring = geom.GetGeometryRef(ring) points = linear_ring.GetPoints() if points: path.moveTo(points[0][0], points[0][1]) for p in points[1:]: path.lineTo(p[0], p[1]) path.closeSubpath() elif geom.GetGeometryType() == ogr.wkbLineString: points = geom.GetPoints() if points: path.moveTo(points[0][0], points[0][1]) for p in points[1:]: path.lineTo(p[0], p[1]) elif geom.GetGeometryType() == ogr.wkbPoint: x, y = geom.GetX(), geom.GetY() path.addEllipse(x - 2, y - 2, 4, 4) return path def auto_zoom(self): """自动缩放视图到数据范围并放大2倍""" if not self.total_bounds.isValid(): return # 设置场景边界 self.scene.setSceneRect(self.total_bounds) # 获取视图可视区域 view_rect = self.graphicsView.viewport().rect() # 计算缩放比例(自动适应 + 2倍放大) transform = QTransform() transform.scale(2, 2) # 先放大2倍 # 应用缩放并居中 self.graphicsView.setTransform(transform) self.graphicsView.fitInView(self.total_bounds, Qt.KeepAspectRatio) # 新增缓冲区分析方法 def buffer_analysis(self): """执行缓冲区分析""" if not hasattr(self, 'current_vector_path'): QMessageBox.warning(self, "警告", "请先打开矢量数据文件!") return # 获取缓冲距离 distance, ok = QInputDialog.getDouble( self, "缓冲区分析", "输入缓冲距离(单位与数据坐标系一致):", 0.0, 0 ) if not ok: return try: # 重新打开数据源获取几何 data_source = ogr.Open(self.current_vector_path, 0) layer = data_source.GetLayer(0) # 创建缓冲区路径 buffer_path = QPainterPath() pen = QPen(Qt.red, 2, Qt.DashLine) brush = QBrush(QColor(255, 0, 0, 50)) # 半透明红色填充 for feature in layer: geom = feature.GetGeometryRef() buffer_geom = geom.Buffer(distance) path = self.geometry_to_qpainterpath(buffer_geom) buffer_path.addPath(path) # 添加到场景 self.scene.addPath(buffer_path, pen, brush) # 更新视图边界 if buffer_path.boundingRect().isValid(): self.total_bounds = self.total_bounds.united(buffer_path.boundingRect()) self.auto_zoom() data_source = None except Exception as e: QMessageBox.critical(self, "错误", f"缓冲区分析失败: {str(e)}") def load_vector_data(self, file_path): self.scene.clear() self.total_bounds = QRectF() try: data_source = ogr.Open(file_path, 0) layer = data_source.GetLayer(0) # 获取字段定义 layer_defn = layer.GetLayerDefn() field_names = [layer_defn.GetFieldDefn(i).GetName() for i in range(layer_defn.GetFieldCount())] # 新增代码:获取并存储坐标系信息 srs = layer.GetSpatialRef() if srs: self.scene.srs = srs # 存储坐标系信息到场景对象 for feature in layer: geom = feature.GetGeometryRef() path = self.geometry_to_qpainterpath(geom) # 创建属性字典 attributes = { "FID": feature.GetFID(), **{name: feature.GetField(name) for name in field_names} } # 使用自定义图形项 item = FeatureItem(path, attributes) item.setPen(QPen(Qt.blue, 1)) item.setBrush(QBrush(Qt.cyan)) item.setZValue(100) # 新增:确保要素在最上层 self.scene.addItem(item) if path.boundingRect().isValid(): self.total_bounds = self.total_bounds.united(path.boundingRect()) data_source = None except Exception as e: print(f"加载失败: {str(e)}") self.current_vector_path = file_path data_source = None # 新增属性查询方法 def enable_query_mode(self): """启用属性查询模式""" self.is_query_mode = not self.is_query_mode self.actionQuery_Attribute.setText("退出查询" if self.is_query_mode else "属性查询") self.graphicsView.setCursor(Qt.CrossCursor if self.is_query_mode else Qt.ArrowCursor) # 新增鼠标事件处理 def mousePressEvent(self, event): # 优先处理空间查询模式 if self.is_spatial_query_mode and event.button() == Qt.LeftButton: scene_pos = self.graphicsView.mapToScene(event.pos()) items = self.scene.items(scene_pos, Qt.IntersectsItemShape, Qt.DescendingOrder) for item in items: if isinstance(item, FeatureItem): # 获取空间关系选择 relations = ["相交", "包含", "被包含", "接触", "重叠"] relation, ok = QInputDialog.getItem( self, "空间关系选择", "请选择空间关系:", relations, 0, False ) if ok: self.perform_spatial_query(item, relation) return return # 处理属性查询模式 if self.is_query_mode and event.button() == Qt.LeftButton: scene_pos = self.graphicsView.mapToScene(event.pos()) items = self.scene.items(scene_pos, Qt.IntersectsItemShape, Qt.DescendingOrder) for item in items: if isinstance(item, FeatureItem): # 确保有属性数据 if hasattr(item, 'attributes') and item.attributes: info = "\n".join([f"{k}: {v}" for k, v in item.attributes.items()]) QMessageBox.information(self, "要素属性", info) else: QMessageBox.warning(self, "提示", "该要素无属性信息") return # 未点击到要素时提示 QMessageBox.warning(self, "提示", "未选中任何要素") return # 其他情况交给父类处理 super().mousePressEvent(event) def draw_centroids(self): """独立质心绘制功能""" if not hasattr(self, 'current_vector_path'): QMessageBox.warning(self, "警告", "请先打开矢量数据文件!") return # 清除已有质心 for item in self.centroid_items: self.scene.removeItem(item) self.centroid_items.clear() try: data_source = ogr.Open(self.current_vector_path, 0) layer = data_source.GetLayer(0) for feature in layer: geom = feature.GetGeometryRef() centroid = geom.Centroid() if centroid: # 获取当前视图变换矩阵 view_transform = self.graphicsView.transform() # 将地理坐标转换为场景坐标 scene_pos = view_transform.map(QPointF( centroid.GetX(), centroid.GetY() )) # 创建质心图形项(使用场景坐标) path = QPainterPath() path.addEllipse( QRectF( scene_pos.x() - 2, # 固定2像素半径 scene_pos.y() - 2, 4, 4 ) ) item = self.scene.addPath( path, QPen(Qt.red, 0.1), QBrush(Qt.red) ) self.centroid_items.append(item) data_source = None self.auto_zoom() except Exception as e: QMessageBox.critical(self, "错误", f"质心绘制失败: {str(e)}") # 新增空间查询模式切换方法 def enable_spatial_query_mode(self): """启用空间查询模式""" self.is_spatial_query_mode = not self.is_spatial_query_mode self.actionSpatial_Query.setText("退出空间查询" if self.is_spatial_query_mode else "空间查询") self.graphicsView.setCursor(Qt.CrossCursor if self.is_spatial_query_mode else Qt.ArrowCursor) if not self.is_spatial_query_mode: self.clear_spatial_query_results() # 新增空间查询核心方法 def perform_spatial_query(self, source_item, relation): """执行空间查询并高亮结果""" self.clear_spatial_query_results() try: # 获取源要素几何 source_geom = self.item_to_geometry(source_item) if not source_geom: return # 获取所有要素 all_items = [item for item in self.scene.items() if isinstance(item, FeatureItem)] # 遍历检查空间关系 for target_item in all_items: target_geom = self.item_to_geometry(target_item) if not target_geom: continue # 执行空间关系判断 if relation == "相交" and source_geom.Intersects(target_geom): self.highlight_item(target_item) elif relation == "包含" and source_geom.Contains(target_geom): self.highlight_item(target_item) elif relation == "被包含" and target_geom.Contains(source_geom): self.highlight_item(target_item) elif relation == "接触" and source_geom.Touches(target_geom): self.highlight_item(target_item) elif relation == "重叠" and source_geom.Overlaps(target_geom): self.highlight_item(target_item) except Exception as e: QMessageBox.critical(self, "错误", f"空间查询失败: {str(e)}") # 新增辅助方法 def item_to_geometry(self, item): """将图形项转换为OGR几何对象""" path = item.path() elements = path.toSubpathPolygons(QTransform()) if not elements: return None # 创建多边形几何 geom = ogr.Geometry(ogr.wkbPolygon) ring = ogr.Geometry(ogr.wkbLinearRing) for point in elements[0]: ring.AddPoint(point.x(), point.y()) ring.CloseRings() geom.AddGeometry(ring) return geom def highlight_item(self, item): """高亮显示查询结果""" original_pen = item.pen() highlight_pen = QPen(Qt.yellow, 3) item.setPen(highlight_pen) self.spatial_query_results.append((item, original_pen)) def clear_spatial_query_results(self): """清除查询结果高亮""" for item, original_pen in self.spatial_query_results: item.setPen(original_pen) self.spatial_query_results.clear() def open_raster_data(self): """打开栅格数据文件""" file_path, _ = QFileDialog.getOpenFileName( self, "打开栅格文件", "", "GeoTIFF (*.tif);;JPEG (*.jpg *.jpeg);;PNG (*.png);;All Files (*)" ) if file_path: try: self.load_raster_data(file_path) self.auto_zoom() except Exception as e: QMessageBox.critical(self, "错误", f"加载栅格失败: {str(e)}") def load_raster_data(self, file_path): """加载栅格数据到视图""" # 打开栅格文件(需要用户修改路径的部分) dataset = gdal.Open(file_path) # 相对路径示例:"./data/raster.tif" # 读取第一个波段 band = dataset.GetRasterBand(1) width = dataset.RasterXSize height = dataset.RasterYSize # 转换为numpy数组 data = band.ReadAsArray() # 创建QImage(注意数据类型转换) if data.dtype == np.uint8: format = QImage.Format.Format_Grayscale8 else: format = QImage.Format.Format_ARGB32 q_img = QImage(data.tobytes(), width, height, format) # 创建像素图项 pixmap = QPixmap.fromImage(q_img) raster_item = self.scene.addPixmap(pixmap) # 处理地理坐标(如果存在) geotransform = dataset.GetGeoTransform() if geotransform: # 计算四个角的坐标 x_origin = geotransform[0] y_origin = geotransform[3] pixel_width = geotransform[1] pixel_height = geotransform[5] # 更新场景边界 x_min = x_origin x_max = x_origin + pixel_width * width y_min = y_origin + pixel_height * height y_max = y_origin self.total_bounds = QRectF( QPointF(x_min, y_min), QPointF(x_max, y_max) ) dataset = None # 关闭数据集 def open_band_combination(self): if not hasattr(self, 'current_raster_path'): QMessageBox.warning(self, "警告", "请先打开栅格数据文件!") return # 复用open_raster_data的逻辑 self.open_raster_data() def open_raster_data(self): file_path, _ = QFileDialog.getOpenFileName( self, "打开栅格文件", "", "GeoTIFF (*.tif);;JPEG (*.jpg *.jpeg);;PNG (*.png);;All Files (*)" ) if file_path: try: dataset = gdal.Open(file_path) num_bands = dataset.RasterCount # 获取用户输入的波段组合 red_band, ok1 = QInputDialog.getInt( self, "波段选择", f"红通道波段号 (1-{num_bands}):", 1, 1, num_bands ) green_band, ok2 = QInputDialog.getInt( self, "波段选择", f"绿通道波段号 (1-{num_bands}):", min(2, num_bands), 1, num_bands ) blue_band, ok3 = QInputDialog.getInt( self, "波段选择", f"蓝通道波段号 (1-{num_bands}):", min(3, num_bands), 1, num_bands ) if not (ok1 and ok2 and ok3): return self.load_raster_data(file_path, red_band, green_band, blue_band) self.auto_zoom() self.current_raster_path = file_path # 新增存储当前路径 except Exception as e: QMessageBox.critical(self, "错误", f"加载栅格失败: {str(e)}") def load_raster_data(self, file_path, red_band=1, green_band=2, blue_band=3): """加载栅格数据到视图(支持波段组合)""" dataset = gdal.Open(file_path) width = dataset.RasterXSize height = dataset.RasterYSize # 读取三个波段数据 def read_band(band_num): band = dataset.GetRasterBand(band_num) data = band.ReadAsArray() # 自动拉伸到0-255范围 data_min = data.min() data_max = data.max() return np.clip(((data - data_min) / (data_max - data_min) * 255), 0, 255).astype(np.uint8) # 合并波段 rgb_array = np.dstack([ read_band(red_band), read_band(green_band), read_band(blue_band) ]) # 创建QImage q_img = QImage( rgb_array.data, width, height, 3 * width, # 每像素3字节(RGB) QImage.Format.Format_RGB888 ) # 创建像素图项 pixmap = QPixmap.fromImage(q_img) self.scene.addPixmap(pixmap) # 处理地理坐标(保持原有逻辑) geotransform = dataset.GetGeoTransform() if geotransform: x_origin = geotransform[0] y_origin = geotransform[3] pixel_width = geotransform[1] pixel_height = geotransform[5] x_min = x_origin x_max = x_origin + pixel_width * width y_min = y_origin + pixel_height * height # 计算下边界 y_max = y_origin # 上边界 # 确保坐标顺序正确 if x_min > x_max: x_min, x_max = x_max, x_min if y_min > y_max: y_min, y_max = y_max, y_min self.total_bounds = QRectF(QPointF(x_min, y_min), QPointF(x_max, y_max)) dataset = None # 新增栅格裁剪方法(必须缩进在类内部) def clip_raster(self): """执行栅格裁剪功能""" if not hasattr(self, 'current_raster_path'): QMessageBox.warning(self, "警告", "请先打开栅格数据文件!") return # 选择裁剪矢量文件 vector_path, _ = QFileDialog.getOpenFileName( self, "选择裁剪区域文件", "", "Shapefile (*.shp);;GeoJSON (*.geojson);;All Files (*)" ) if not vector_path: return try: # 获取原始栅格信息 src_ds = gdal.Open(self.current_raster_path) geotransform = src_ds.GetGeoTransform() proj = src_ds.GetProjection() # 获取矢量范围 vector_ds = ogr.Open(vector_path) layer = vector_ds.GetLayer() feature = layer.GetNextFeature() geom = feature.GetGeometryRef() x_min, x_max, y_min, y_max = geom.GetEnvelope() # 创建临时裁剪结果文件 import os # 确保导入os模块 output_path = os.path.splitext(self.current_raster_path)[0] + "_clipped.tif" # 执行裁剪操作 options = gdal.WarpOptions( format='GTiff', outputBounds=[x_min, y_min, x_max, y_max], dstSRS=proj ) gdal.Warp(output_path, src_ds, options=options) # 加载裁剪结果 self.load_raster_data(output_path) self.auto_zoom() # 清理资源 src_ds = None vector_ds = None except Exception as e: QMessageBox.critical(self, "错误", f"栅格裁剪失败: {str(e)}") # 新增波段运算方法 def band_calculation(self): """执行波段运算(示例为NDVI计算)""" if not hasattr(self, 'current_raster_path'): QMessageBox.warning(self, "警告", "请先打开栅格数据文件!") return try: # 获取用户输入参数 red_band, ok1 = QInputDialog.getInt( self, "波段选择", "输入红波段编号 (1-based):", 1, 1, 100 ) nir_band, ok2 = QInputDialog.getInt( self, "波段选择", "输入近红外波段编号 (1-based):", 4, 1, 100 ) if not (ok1 and ok2): return # 读取栅格数据 dataset = gdal.Open(self.current_raster_path) red = dataset.GetRasterBand(red_band).ReadAsArray() nir = dataset.GetRasterBand(nir_band).ReadAsArray() # 执行NDVI计算 ndvi = np.where( (nir + red) == 0, 0, (nir - red) / (nir + red) ).astype(np.float32) # 创建输出文件 output_path, _ = QFileDialog.getSaveFileName( self, "保存结果", "", "GeoTIFF (*.tif)" ) if not output_path: return # 写入结果 driver = gdal.GetDriverByName('GTiff') out_ds = driver.Create( output_path, dataset.RasterXSize, dataset.RasterYSize, 1, gdal.GDT_Float32 ) out_ds.SetGeoTransform(dataset.GetGeoTransform()) out_ds.SetProjection(dataset.GetProjection()) out_ds.GetRasterBand(1).WriteArray(ndvi) out_ds.FlushCache() # 清理资源 dataset = None out_ds = None QMessageBox.information(self, "成功", f"NDVI计算结果已保存至:\n{output_path}") except Exception as e: QMessageBox.critical(self, "错误", f"波段运算失败: {str(e)}") # 新增保存方法 def save_vector_data(self): """保存当前场景中的矢量数据""" if not hasattr(self, 'current_vector_path'): QMessageBox.warning(self, "警告", "没有可保存的矢量数据!") return # 获取保存路径 file_path, _ = QFileDialog.getSaveFileName( self, "保存矢量文件", "", "Shapefile (*.shp);;GeoJSON (*.geojson);;All Files (*)" ) if not file_path: return try: # 获取原始数据源信息 src_ds = ogr.Open(self.current_vector_path) src_layer = src_ds.GetLayer(0) src_defn = src_layer.GetLayerDefn() # 创建目标数据源 driver = ogr.GetDriverByName("ESRI Shapefile" if file_path.endswith(".shp") else "GeoJSON") if os.path.exists(file_path): driver.DeleteDataSource(file_path) dst_ds = driver.CreateDataSource(file_path) # 创建图层(保持与原始数据相同的坐标系) dst_layer = dst_ds.CreateLayer( "saved_features", srs=src_layer.GetSpatialRef(), geom_type=ogr.wkbPolygon ) # 复制字段定义 for i in range(src_defn.GetFieldCount()): field_defn = src_defn.GetFieldDefn(i) dst_layer.CreateField(field_defn) # 遍历场景中的所有要素项 for item in self.scene.items(): if isinstance(item, FeatureItem): # 创建新要素 feature = ogr.Feature(dst_layer.GetLayerDefn()) # 复制属性 for key, value in item.attributes.items(): if key == "FID": continue # FID通常自动生成 if feature.GetFieldIndex(key) != -1: feature.SetField(key, str(value)) # 转换几何 geom = self.item_to_geometry(item) if geom: feature.SetGeometry(geom) dst_layer.CreateFeature(feature) feature = None # 添加缓冲区要素(如果存在) self.save_additional_features(dst_layer, "buffer") # 添加质心要素(如果存在) self.save_additional_features(dst_layer, "centroid") dst_ds = None src_ds = None QMessageBox.information(self, "成功", f"数据已保存至:\n{file_path}") except Exception as e: QMessageBox.critical(self, "错误", f"保存失败: {str(e)}") # 新增辅助保存方法 def save_additional_features(self, layer, feature_type): """保存附加要素(缓冲区/质心)""" items = [] if feature_type == "buffer": items = [item for item in self.scene.items() if item.pen().style() == Qt.DashLine and item.pen().color() == Qt.red] elif feature_type == "centroid": items = self.centroid_items for item in items: geom = self.item_to_geometry(item) if geom: feature = ogr.Feature(layer.GetLayerDefn()) feature.SetGeometry(geom) # 添加类型标识字段 feature.SetField("FEATURE_TYPE", feature_type.upper()) layer.CreateFeature(feature) feature = None # 在item_to_geometry方法中添加点要素支持 def item_to_geometry(self, item): """增强版几何转换(支持点要素)""" path = item.path() elements = path.toSubpathPolygons(QTransform()) if not elements: # 处理点要素 if isinstance(item, QGraphicsPathItem): path = item.path() if path.elementCount() == 1 and path.elementAt(0).isMoveTo(): pt = path.elementAt(0) geom = ogr.Geometry(ogr.wkbPoint) geom.AddPoint(pt.x, pt.y) return geom return None # 原有多边形处理逻辑 geom = ogr.Geometry(ogr.wkbPolygon) ring = ogr.Geometry(ogr.wkbLinearRing) for point in elements[0]: ring.AddPoint(point.x(), point.y()) ring.CloseRings() geom.AddGeometry(ring) return geom # 新增栅格保存方法 def save_raster_data(self): """保存当前显示的栅格数据""" if not hasattr(self, 'current_raster_path'): QMessageBox.warning(self, "警告", "没有可保存的栅格数据!") return # 获取保存路径和格式 file_path, selected_filter = QFileDialog.getSaveFileName( self, "保存栅格文件", "", "GeoTIFF (*.tif);;JPEG (*.jpg *.jpeg);;PNG (*.png);;All Files (*)" ) if not file_path: return try: # 获取当前显示的栅格数据 dataset = gdal.Open(self.current_raster_path) band_count = dataset.RasterCount width = dataset.RasterXSize height = dataset.RasterYSize # 创建输出数据集 driver = gdal.GetDriverByName(self.get_driver_name(file_path)) out_ds = driver.Create( file_path, width, height, band_count, self.get_gdal_datatype(dataset) ) # 复制地理参考信息 out_ds.SetGeoTransform(dataset.GetGeoTransform()) out_ds.SetProjection(dataset.GetProjection()) # 复制所有波段数据 for band_num in range(1, band_count + 1): in_band = dataset.GetRasterBand(band_num) out_band = out_ds.GetRasterBand(band_num) out_band.WriteArray(in_band.ReadAsArray()) out_band.FlushCache() # 清理资源 out_ds = None dataset = None QMessageBox.information(self, "成功", f"栅格数据已保存至:\n{file_path}") except Exception as e: QMessageBox.critical(self, "错误", f"栅格保存失败: {str(e)}") # 新增辅助方法 def get_driver_name(self, file_path): """根据文件扩展名获取GDAL驱动名称""" ext = os.path.splitext(file_path)[1].lower() return { '.tif': 'GTiff', '.jpg': 'JPEG', '.jpeg': 'JPEG', '.png': 'PNG' }.get(ext, 'GTiff') def get_gdal_datatype(self, dataset): """获取与原始数据集匹配的GDAL数据类型""" sample_band = dataset.GetRasterBand(1) return { gdal.GDT_Byte: gdal.GDT_Byte, gdal.GDT_UInt16: gdal.GDT_UInt16, gdal.GDT_Int16: gdal.GDT_Int16, gdal.GDT_UInt32: gdal.GDT_UInt32, gdal.GDT_Int32: gdal.GDT_Int32, gdal.GDT_Float32: gdal.GDT_Float32, gdal.GDT_Float64: gdal.GDT_Float64 }.get(sample_band.DataType, gdal.GDT_Float32)
最新发布
07-18
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值