arcpy批量重分类,分区统计

批量重分类,先去gis里记住分类间断点。

import arcpy
import os
from arcpy.sa import *
import time

# 检查空间分析扩展许可
arcpy.CheckOutExtension("Spatial")

# 设置工作环境
input_path = r'D:\Aganhan\AYN\云南综合风险栅格'
output_path = r'D:\Adriving\定量分析\soc_reclass'

# 记录开始时间
start_time = time.time()

# 如果输出文件夹不存在则创建
if not os.path.exists(output_path):
    os.makedirs(output_path)
    print(f"创建输出文件夹: {output_path}")

# 定义重分类规则
reclass_range = []
reclass_range.append([0, 0.02, 1])           # 低风险
reclass_range.append([0.02, 0.05, 2])        # 较低风险
reclass_range.append([0.05, 0.1, 3])         # 中风险
reclass_range.append([0.1, 0.2, 4])          # 较高风险
reclass_range.append([0.2, 9999, 5])         # 高风险

# 创建重分类对象
remap = arcpy.sa.RemapRange(reclass_range)

# 获取输入文件夹中的所有tif文件
tif_files = [f for f in os.listdir(input_path) if f.endswith('.tif')]
total_files = len(tif_files)

print(f"\n开始处理,共有 {total_files} 个文件需要处理")
print("----------------------------------------")

# 处理每个tif文件
for index, tif in enumerate(tif_files, 1):
    try:
        # 获取年份
        year = tif.split('.')[0][-4:]
        
        # 构建输入和输出文件路径
        input_raster = os.path.join(input_path, tif)
        output_raster = os.path.join(output_path, f'风险等级_{year}.tif')
        
        # 计算处理进度
        progress = (index / total_files) * 100
        
        print(f"\n处理第 {index}/{total_files} 个文件 ({progress:.1f}%)")
        print(f"年份: {year}")
        print(f"输入文件: {tif}")
        print(f"输出文件: 风险等级_{year}.tif")
        
        # 记录单个文件处理开始时间
        file_start_time = time.time()
        
        # 执行重分类
        print("执行重分类...")
        outReclass = Reclassify(input_raster, "VALUE", remap)
        
        # 保存结果
        print("保存结果...")
        outReclass.save(output_raster)
        
        # 计算单个文件处理时间
        file_time = time.time() - file_start_time
        print(f"完成处理")
        print(f"处理时间: {file_time:.2f} 秒")
        
    except Exception as e:
        print(f"处理 {tif} 时出错: {str(e)}")

# 计算总处理时间
total_time = time.time() - start_time
hours = int(total_time // 3600)
minutes = int((total_time % 3600) // 60)
seconds = total_time % 60

print("\n----------------------------------------")
print("所有文件处理完成!")
print(f"总处理时间: {hours}小时 {minutes}分钟 {seconds:.2f}秒")
print(f"处理结果保存在: {output_path}")

# 释放许可
arcpy.CheckInExtension("Spatial")

分区统计工具,定量分析影响的人口和经济:

import arcpy
import os
import pandas as pd
from arcpy.sa import *
import time

# 检查空间分析扩展许可
arcpy.CheckOutExtension("Spatial")

# 设置工作环境
reclass_path = r'D:\Adriving\定量分析\soc_reclass'
gdp_path = r'D:\Adriving\GDP'
pop_path = r'D:\Adriving\人口数据'
output_path = r'D:\Adriving\定量分析\统计结果'

# 如果输出文件夹不存在则创建
if not os.path.exists(output_path):
    os.makedirs(output_path)
    print(f"创建输出文件夹: {output_path}")

# 创建存储结果的列表
results = []

# 获取所有重分类后的tif文件
reclass_files = [f for f in os.listdir(reclass_path) if f.endswith('.tif')]
total_files = len(reclass_files)

print(f"\n开始处理,共有 {total_files} 个年份需要处理")
print("----------------------------------------")

# 处理每个年份
for index, reclass_file in enumerate(reclass_files, 1):
    try:
        # 获取年份
        year = reclass_file.split('_')[1].split('.')[0]
        
        # 构建文件路径
        reclass_raster = os.path.join(reclass_path, reclass_file)
        gdp_raster = os.path.join(gdp_path, f'【立方数据学社】云南省{year}GDP.tif')
        pop_raster = os.path.join(pop_path, f'学术严选_china_{year}.tif')
        
        print(f"\n处理第 {index}/{total_files} 个年份: {year}")
        
        # 确保重分类栅格为整型
        reclass_int = Int(Raster(reclass_raster))
        
        # 使用分区统计工具
        print(f"统计 {year} 年数据...")
        
        # 统计GDP
        print("计算GDP...")
        gdp_stats = ZonalStatisticsAsTable(reclass_int, "Value", gdp_raster, 
                                         os.path.join(output_path, f"gdp_stats_{year}.dbf"), 
                                         "DATA", "SUM")
        
        # 统计人口
        print("计算人口...")
        pop_stats = ZonalStatisticsAsTable(reclass_int, "Value", pop_raster, 
                                         os.path.join(output_path, f"pop_stats_{year}.dbf"), 
                                         "DATA", "SUM")
        
        # 读取统计结果
        with arcpy.da.SearchCursor(gdp_stats, ["Value", "SUM"]) as gdp_cursor:
            for gdp_row in gdp_cursor:
                risk_level = gdp_row[0]
                gdp_sum = gdp_row[1]
                
                # 获取对应的人口数据
                pop_sum = 0
                with arcpy.da.SearchCursor(pop_stats, ["Value", "SUM"]) as pop_cursor:
                    for pop_row in pop_cursor:
                        if pop_row[0] == risk_level:
                            pop_sum = pop_row[1]
                            break
                
                results.append({
                    'Year': year,
                    'Risk_Level': risk_level,
                    'GDP_Sum': gdp_sum,
                    'Population_Sum': pop_sum
                })
        
        print(f"完成 {year} 年统计")
        
    except Exception as e:
        print(f"处理 {year} 年数据时出错: {str(e)}")
        import traceback
        print(traceback.format_exc())

# 将结果转换为DataFrame并保存
print("\n保存统计结果...")
df = pd.DataFrame(results)
excel_file = os.path.join(output_path, '风险分区统计结果.xlsx')
df.to_excel(excel_file, index=False)

print("\n----------------------------------------")
print("所有统计完成!")
print(f"统计结果已保存至: {excel_file}")

# 释放许可
arcpy.CheckInExtension("Spatial")

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值