# -*- coding: utf-8 -*- __author__ = 'wanger' __description__ = '批量生成栅格数据元数据相关工具' __date__ = '2024-11-25' __copyright__ = '(C) 2024 by siwei' __revision__ = '1.0' import os from osgeo import osr from qgis._core import QgsCoordinateReferenceSystem from qgis.core import ( QgsRasterLayer, QgsVectorLayer, QgsSpatialIndex, QgsFeatureRequest, QgsGeometry ) # 通过文件路径来获取文件名和文件后缀 def getImageNumberAndFormat(filepath): # 获取文件名(包括后缀) file_name = os.path.basename(filepath) # 获取文件名(不包括后缀) file_name_no_ext = os.path.splitext(file_name)[0] # 获取文件后缀 file_ext = os.path.splitext(file_name)[1][1:] return file_name_no_ext, file_ext # 通过文件路径获取文件大小 def getFileSize(filepath): try: file_size = os.path.getsize(filepath) size_mb = file_size / (1024 * 1024) return f"{size_mb:.2f}" except FileNotFoundError: print("文件未找到,请检查路径是否正确!") return 0 def getImageLeftTop(dataset): # 获取地理变换参数 geotransform = dataset.GetGeoTransform() if not geotransform: raise ValueError("无法获取地理变换参数!") # 左上角像素的坐标(投影坐标系) x_min = geotransform[0] y_max = geotransform[3] # 获取栅格的空间参考(投影坐标系) spatial_ref = osr.SpatialReference() spatial_ref.ImportFromWkt(dataset.GetProjection()) # 转换为地理坐标系(经纬度) if spatial_ref.IsProjected(): # spatial_ref_geo = spatial_ref.CloneGeogCS() wgs84_crs = osr.SpatialReference() wgs84_crs.ImportFromEPSG(4326) transform = osr.CoordinateTransformation(spatial_ref, wgs84_crs) lon, lat, _ = transform.TransformPoint(x_min, y_max) else: lon, lat = x_min, y_max # 如果已经是地理坐标系,直接使用 return { "lon": lat, "lat": lon, "y": int(x_min), "x": int(y_max) } # 获取中央子午线 def getCentralMeridian(spatial_ref): proj4 = spatial_ref.ExportToProj4() print(f"Projection PROJ.4: {proj4}") central_meridian = None for param in proj4.split(): if param.startswith('+lon_0'): central_meridian = param.split('=')[1] break if central_meridian: return central_meridian else: print("No central meridian found.") return None # 通过栅格数据与矢量数据叠加分析 获取相交面积最大的一块 返回要素 def get_largest_overlapping_feature(raster_path, shp_path): """ 查询与栅格相交的矢量要素中面积最大的一块,并返回指定字段的值。 :param raster_path: 栅格文件路径 :param shp_path: 矢量文件路径 :return: 面积最大的要素 """ # 加载栅格图层 raster_layer = QgsRasterLayer(raster_path, "Raster Layer") if not raster_layer.isValid(): raise Exception(f"栅格文件无效: {raster_path}") # 加载矢量图层 vector_layer = QgsVectorLayer(shp_path, "Vector Layer", "ogr") if not vector_layer.isValid(): raise Exception(f"SHP文件无效: {shp_path}") # 获取栅格的范围 raster_extent = raster_layer.extent() # 创建空间索引以加速查询 spatial_index = QgsSpatialIndex(vector_layer.getFeatures()) # 查询与栅格范围相交的要素ID intersecting_ids = spatial_index.intersects(raster_extent) # 遍历这些要素并计算相交区域的面积 max_area = 0 largest_feature = None for fid in intersecting_ids: feature = vector_layer.getFeature(fid) intersection = feature.geometry().intersection(QgsGeometry.fromRect(raster_extent)) # 检查是否有有效相交区域 if intersection.isEmpty(): continue # 计算相交区域的面积 area = intersection.area() # 找到面积最大的要素 if area > max_area: max_area = area largest_feature = feature if largest_feature: return largest_feature else: return None # 处理数值 def processingNumericalValues(value): print(value) if isinstance(value, str): try: value = float(value) except ValueError: print(f"无法将 '{value}' 转换为浮动数") return None if value.is_integer(): return int(value) else: return round(value, 2)