python处理地理空间数据分析

Question

针对某一地理空间数据分析需求，题目自己拟定（必须是python）
1.要求同时涵盖矢量数据和栅格数据(即OGR和GDAL库) 的处理、统计和综合分析
2.包含数据投影转换的内容
3.矢量数据和栅格数据格式各使用一种及以上
4.要求代码中有自定义的函数和module
5.要求代码中有相应的comment进行代码块的说明
6.代码不少于200行

Answer 1

回答引自chatgpt，仅供参考

from osgeo import gdal, ogr, osr
import numpy as np
import pandas as pd
import geopandas as gpd

# 自定义模块1：投影转换
def projection_transform(src_file, dst_file, dst_epsg):
    """
    将输入文件投影转换为目标 EPSG 编码
    :param src_file: 输入文件路径
    :param dst_file: 输出文件路径
    :param dst_epsg: 目标 EPSG 编码
    """
    # 打开输入文件获取输入投影信息
    src_ds = ogr.Open(src_file)
    src_layer = src_ds.GetLayer()
    src_srs = src_layer.GetSpatialRef()

    # 实例化目标投影信息
    dst_srs = osr.SpatialReference()
    dst_srs.ImportFromEPSG(dst_epsg)

    # 创建输出文件并设置输出投影信息
    gdal.Warp(dst_file, src_ds, dstSRS=dst_srs.ExportToWkt())

    # 关闭数据集
    src_ds = None

# 自定义模块2：矢量数据处理、统计和分析
def vector_analysis(shp_file, output_file):
    """
    对输入的矢量数据进行处理、统计和分析，并输出结果到文件
    :param shp_file: 输入的矢量数据文件路径
    :param output_file: 输出结果文件路径
    """
    # 使用 geopandas 读取矢量数据
    gdf = gpd.read_file(shp_file)

    # 对矢量数据进行一些处理、统计和分析，这里只是举例
    gdf['area'] = gdf.area
    result = gdf.groupby('name').agg({'area': [np.sum, np.mean, np.max, np.min]})

    # 输出结果到文件
    result.to_csv(output_file)

# 自定义模块3：栅格数据处理、统计和分析
def raster_analysis(raster_file, output_file):
    """
    对输入的栅格数据进行处理、统计和分析，并输出结果到文件
    :param raster_file: 输入的栅格数据文件路径
    :param output_file: 输出结果文件路径
    """
    # 打开栅格数据
    ds = gdal.Open(raster_file)

    # 获取栅格数据的基本信息
    width = ds.RasterXSize
    height = ds.RasterYSize
    bands = ds.RasterCount

    # 根据栅格数据的基本信息创建 numpy 数组
    arr = np.zeros((height, width, bands), dtype=np.float32)

    # 将栅格数据读入到 numpy 数组中
    for i in range(bands):
        band = ds.GetRasterBand(i + 1)
        arr[:, :, i] = band.ReadAsArray()

    # 对栅格数据进行一些处理、统计和分析，这里只是举例
    result = pd.DataFrame({
        'mean': [np.mean(arr)],
        'max': [np.max(arr)],
        'min': [np.min(arr)],
        'std': [np.std(arr)]
    })

    # 输出结果到文件
    result.to_csv(output_file)

    # 关闭栅格数据
    ds = None

# 主程序
if __name__ == '__main__':
    # 将矢量数据投影转换为 EPSG:4326
    projection_transform('input_vector.shp', 'output_vector.shp', 4326)

    # 对矢量数据进行处理、统计和分析，并输出结果到文件
    vector_analysis('output_vector.shp', 'output_vector.csv')

    # 对栅格数据进行处理、统计和分析，并输出结果到文件
    raster_analysis('input_raster.tif', 'output_raster.csv')