Python地理数据处理 27:基于Arcpy批量处理已矫正的worldclim2.1未来气候数据——投影、重采样、多波段拆分以及裁剪

news2024/10/1 20:57:32

Arcpy批量处理已矫正的worldclim2.1未来气候数据

  • 1. 写在前面
  • 2.实现代码

1. 写在前面

  前面我写了一篇关于如何使用ArcGIS自带的Python工具处理worldclim数据的多波段数据的文章,而这只是处理该数据的其中一步。要想得到满足要求的数据,还需要其他操作,依次为投影为指定投影坐标系(Albers)、重采样为1000m空间分辨率、多波段拆分以及裁剪。
  今天我把以上所有的相关操作基于一套代码进行演示,可以实现一键操作。并且我使用的是已矫正的worldclim2.1未来气候数据,该数据排除了异常值。

2.实现代码

  使用以上代码之前,你需要建立一个文件夹,例如:C:\Users\Jacksonzhao\Desktop\sspdata\ssp126。再将worldclim数据放入其中,这里我选择了worldclim的ssp126_2030s和ssp126_2050s两幅影像数据进行处理,之后便可以运行代码,喝一杯茶的功夫就好了:

在这里插入图片描述

# -*- coding: cp936 -*-
import arcpy
import shutil
from arcpy.sa import *
import os
import time

# 检查并获取 Spatial Analyst 扩展许可
arcpy.CheckOutExtension('Spatial')
start_time = time.time()  # 开始时间

# *****************************************************************************************
# 基础路径和子文件夹
subfolder = 'ssp126'  # 需要修改
base_path = os.path.join(r'C:\Users\Jacksonzhao\Desktop\sspdata', subfolder)

# 设置工作空间
arcpy.env.workspace = base_path

# 获取工作空间中的所有tif格式文件
raster_lists = arcpy.ListRasters("*.tif")
print(raster_lists)

for raster_list in raster_lists:
     #print(raster_list)
     raster_list = [raster_list]
     # 检查是否找到了tif文件
     if raster_list:
         # 遍历文件列表
         for raster in raster_list:
             print("*****************************************************************************")
             print("正在处理文件:{}************************************************".format(raster))
             # 移除.tif扩展名
             raster_name = os.path.splitext(raster)[0]
             #print(raster_name)  # 打印不含扩展名的文件名

             # 子文件夹模板列表
             subfolders = [
                 '{}/Albers998m'.format(raster_name),
                 '{}/Albers1000m'.format(raster_name),
                 '{}/Albers1000m_19bio'.format(raster_name),
                 '{}_Finish'.format(raster_name)
             ]
             
             print("=========================文件夹创建=========================")
             # 遍历子文件夹模板列表,为每个raster_name创建子文件夹
             for subfolder_template in subfolders:
                 subfolder_path = subfolder_template.format(raster_name=raster_name)
                 folder_path = os.path.join(base_path, subfolder_path)

                 # 检查文件夹是否已经存在
                 if not os.path.exists(folder_path):
                     
                     os.makedirs(folder_path)
                     print("文件夹 '{}' 已创建".format(os.path.basename(folder_path)))
                 else:
                     print("文件夹 '{}' 已存在".format(os.path.basename(folder_path)))

             print("所有文件夹创建完毕,耗时:{:.2f}分====================================".format((time.time()  - start_time)/60))


             # 1、投影为albers
             print("=========================投影为Albers=========================")
             arcpy.env.workspace = base_path
             #rasters = arcpy.ListRasters("*","tif")
             for raster in raster_list:
                 out = os.path.join(base_path, raster_name, 'Albers998m', raster)  # 输出路径
                 print(os.path.join(base_path, raster_name, 'Albers998m', raster))
                 arcpy.ProjectRaster_management(raster, out, "PROJCS['Asia_North_Albers_WGS84_LCR',GEOGCS['GCS_WGS_1984',DATUM['D_WGS_1984',SPHEROID['WGS_1984',6378137.0,298.257223563]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]],PROJECTION['Albers'],PARAMETER['False_Easting',0.0],PARAMETER['False_Northing',0.0],PARAMETER['Central_Meridian',105.0],PARAMETER['Standard_Parallel_1',25.0],PARAMETER['Standard_Parallel_2',47.0],PARAMETER['Latitude_Of_Origin',0.0],UNIT['Meter',1.0]]")
             print("操作完毕,耗时:{:.2f}分====================================".format((time.time()  - start_time)/60))


             # 2、重采样为1000m
             print("=========================重采样为1000m=========================")
             arcpy.env.workspace = os.path.join(base_path, raster_name, 'Albers998m')

             input_folder = os.path.join(base_path, raster_name, 'Albers998m')
             output_folder = os.path.join(base_path, raster_name, 'Albers1000m')

             resampling_resolution = "1000" # 设置重采样分辨率

             tif_list = arcpy.ListFiles("*.tif")

             # 循环处理每个tif文件
             for tif_file in tif_list:
                 #构建输入和输出路径
                 input_path = os.path.join(input_folder, tif_file)
                 output_path = os.path.join(output_folder, tif_file)
                 # 执行重采样
                 arcpy.Resample_management(input_path, output_path, resampling_resolution, "NEAREST")
                 
                 print(output_path)
             print("操作完毕,耗时:{:.2f}分====================================".format((time.time()  - start_time)/60))



             # 3、拆分为单波段
             print("=========================多波段拆分为单波段=========================")
             arcpy.env.workspace = os.path.join(base_path, raster_name, 'Albers1000m')
             output_directory = os.path.join(base_path, raster_name, 'Albers1000m_19bio')
             print(output_directory)

             # 获取栅格列表
             raster_list = arcpy.ListRasters("*.tif")

             # 遍历栅格集合
             for raster in raster_list:
                 # 获取栅格数据集的全路径
                 raster_full_path = os.path.join(arcpy.env.workspace, raster)
                 
                 # 创建栅格对象
                 raster_obj = arcpy.Raster(raster_full_path)
                 
                 # 获取波段数量
                 number_of_bands = raster_obj.bandCount
                 
                 # 构建基本的文件名
                 base_filename = "wc2.1BCC_" + raster_name
                 
                 # 遍历所有波段
                 for i in range(1, number_of_bands + 1):
                     # 提取单波段并创建栅格图层
                     band_name = "band_" + str(i)
                     arcpy.MakeRasterLayer_management(raster_full_path, band_name, "", "", i)

                     # 保存单波段影像
                     output_path = os.path.join(output_directory, base_filename + "_Bio{0}.tif".format(i))
                     arcpy.CopyRaster_management(band_name, output_path)
                     
                     print("Saved:", output_path)
             print("操作完毕,耗时:{:.2f}分====================================".format(time.time()  - start_time))



             # 4、裁剪
             print("=========================影像裁剪=========================")
             finish_folder = os.path.join(base_path, '{}_Finish'.format(raster_name))
             arcpy.env.workspace = os.path.join(base_path, raster_name, 'Albers1000m_19bio')
             rasters = arcpy.ListRasters("*", "tif")
             mask ="C:/Users/Jacksonzhao/Desktop/CSCD/Data/Bio_Topo/envi1km/wc2.1_30s_bio_1.tif"

             for raster in rasters:
                 print(raster)
                 out = os.path.join(finish_folder,  raster)
                 arcpy.gp.ExtractByMask_sa(raster, mask, out)
                 print("Clip_" + raster + "has done!")
             print("操作完毕,耗时:{:.2f}分====================================".format((time.time()  - start_time)/60))



             # 5、删除文件夹
             folder_paths = [
                  os.path.join(base_path, raster_name, 'Albers998m'),
                  os.path.join(base_path, raster_name, 'Albers1000m_19bio')
                  ]

             # 检查并删除文件夹
             for folder_path in folder_paths:
                  if os.path.isdir(folder_path):
                       shutil.rmtree(folder_path)
                       print("文件夹 : {} 及其内容已删除".format(folder_path))
                  else:
                       print("文件夹 : {} 不存在".format(folder_path))
             end_time = time.time()  # 开始时间
             print("所有操作完毕,耗时:{:.2f}分*****************************************".format((end_time - start_time)/60))

     else:
         print("没有找到tif文件。")

print("所有操作完毕!!!,耗时:{:.2f}分*****************************************".format((end_time - start_time)/60))

处理结果:
在这里插入图片描述

BCC_ssp126_2030s_Bio1结果如下:
在这里插入图片描述

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.coloradmin.cn/o/2183136.html

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈,一经查实,立即删除!

相关文章

自闭症寄宿学校 vs. 日常教育:为孩子提供更多可能

在探索自闭症儿童的教育路径时,家长们往往面临一个重大的选择:是选择传统的日常教育环境,还是寻找专为自闭症儿童设计的寄宿学校?广州的星贝育园自闭症儿童寄宿制学校,以其独特的教育模式和全方位的关怀体系&#xff0…

大数据毕业设计选题推荐-个性化图书推荐系统-Python数据可视化-Hive-Hadoop-Spark

✨作者主页:IT毕设梦工厂✨ 个人简介:曾从事计算机专业培训教学,擅长Java、Python、PHP、.NET、Node.js、GO、微信小程序、安卓Android等项目实战。接项目定制开发、代码讲解、答辩教学、文档编写、降重等。 ☑文末获取源码☑ 精彩专栏推荐⬇…

python中的find函数怎么用

Python find() 方法检测字符串中是否包含子字符串 str ,如果指定 beg(开始) 和 end(结束) 范围,则检查是否包含在指定范围内,如果包含子字符串返回开始的索引值,否则返回-1。 语法 …

Netty系列-6 Netty消息处理流程

背景 前文介绍了Netty服务端的启动流程,服务端启动后可以处理客户端发送的请求,包括连接请求和普通消息。 1.处理连接 当客户端有连接请求到达时,服务器会创建通道并将通道注册到选择器上,处理逻辑与NIO中实现完全一致。 详细流…

虚拟机、ubantu不能连接网络,解决办法

虚拟机、ubantu不能连接网络,解决办法 物理机OS: [Windows10 专业版](https://so.csdn.net/so/search?qWindows10 专业版&spm1001.2101.3001.7020) 虚拟机平台: VMware Workstation 16 Pro 虚拟机OS: Ubuntu 18.04 自动配…

英语音标与重弱读

英语中,比较重要的是音标。但事实上,我们对音标的学习还是比较少的,对它的理解也是比较少的。 一、音标 2个半元音 [w][j] 5个长元音:[i:] [ə:] [ɔ:] [u:] [ɑ:] 7个短元音:[i] [ə] [ɔ] [u] [] [e] [ʌ] 8个双元音…

车辆重识别(2020NIPS去噪扩散概率模型)论文阅读2024/9/27

[2] Denoising Diffusion Probabilistic Models 作者:Jonathan Ho Ajay Jain Pieter Abbeel 单位:加州大学伯克利分校 摘要: 我们提出了高质量的图像合成结果使用扩散概率模型,一类潜变量模型从非平衡热力学的考虑启发。我们的最…

【mmengine】配置器(config)(入门)读取与使用

一、 介绍 MMEngine 实现了抽象的配置类(Config),为用户提供统一的配置访问接口。 配置类能够支持不同格式的配置文件,包括 python,json,yaml,用户可以根据需求选择自己偏好的格式。 配置类提供…

leetcode力扣刷题系列——【座位预约管理系统】

题目 请你设计一个管理 n 个座位预约的系统,座位编号从 1 到 n 。 请你实现 SeatManager 类: SeatManager(int n) 初始化一个 SeatManager 对象,它管理从 1 到 n 编号的 n 个座位。所有座位初始都是可预约的。 int reserve() 返回可以预约座…

单调队列应用介绍

单调队列应用介绍 定义应用场景实现模板具体示例滑动窗口最大值问题描述问题分析代码实现带限制的子序列和问题描述问题分析代码实现跳跃游戏问题描述问题分析代码实现定义 队列(Queue)是另一种操作受限的线性表,只允许元素从队列的一端进,另一端出,具有先进先出(FIFO)的特…

系统信息规划-系统架构师(七十四)

1前驱图 解析: 当S1执行完,C1S2并行执行,C1和S2执行完,P1,C2,S3并行执行,同理,P2C3并行执行。 直接制约则表示C1和P1受S1制约。 间接则代表S2和S3受S1制约。 2系统移植也是系统构建的一种实现方…

学习记录:js算法(五十一):统计二叉树中好节点的数目

文章目录 统计二叉树中好节点的数目网上思路 总结 统计二叉树中好节点的数目 给你一棵根为 root 的二叉树,请你返回二叉树中好节点的数目。 「好节点」X 定义为:从根到该节点 X 所经过的节点中,没有任何节点的值大于 X 的值。 图一&#xff1…

长江存储致态TiPlus7100 4TB满盘读写测试:性能几乎没有下降

一、前言:看看满盘状态下致态TiPlus7100 4TB性能会如何! 现在还有很多同学对于长江存储品牌的存储产品不太信任,在选择SSD时会优先考虑三星、西数这样的品牌。 有鉴于此,我们此次会将手上的长江存储致态TiPlus7100 4TB SSD进行更严…

【STM32单片机_(HAL库)】4-2-1【定时器TIM】定时器输出PWM实现呼吸灯实验

1.硬件 STM32单片机最小系统LED灯模块 2.软件 pwm驱动文件添加定时器HAL驱动层文件添加GPIO常用函数定时器输出PWM配置步骤main.c程序 #include "sys.h" #include "delay.h" #include "led.h" #include "pwm.h"int main(void) {HA…

音视频入门基础:FLV专题(10)——Script Tag实例分析

一、引言 在《音视频入门基础:FLV专题(9)——Script Tag简介》中对FLV文件的Script Tag进行了简介。下面用一个具体的例子来对Script Tag进行分析。 二、Script Tag的Tag header实例分析 用notepad打开《音视频入门基础:FLV专题…

鸿蒙跨端实践-JS虚拟机架构实现

作者:京东科技 杜强强 前言 在Roma跨端方案中,JS虚拟机是框架的核心,负责执行动态化的JS代码。在Android平台采用了基于V8的J2V8,iOS平台则使用了系统自带的JSCore,而在HarmonyOS中,由于业界无类似的框架&a…

C++11_左值引用与右值引用

在C11之前,是没有右值引用的概念的,在C11之后才新增了右值引用。其实无论是左值引用还是右值引用都是给对象取别名。 认识左值和右值 什么是左值? 左值是一个表示数据的表达式(如变量名或解引用的指针),我们可以获取它的地址可…

YOLOv11改进策略【损失函数篇】| Shape-IoU:考虑边界框形状和尺度的更精确度量

一、本文介绍 本文记录的是改进YOLOv11的损失函数,将其替换成Shape-IoU。现有边界框回归方法通常考虑真实GT(Ground Truth)框与预测框之间的几何关系,通过边界框的相对位置和形状计算损失,但忽略了边界框本身的形状和…

PV大题--专题突破

写在前面: PV大题考查使用伪代码控制进程之间的同步互斥关系,它需要我们一定的代码分析能力,算法设计能力,有时候会给你一段伪代码让你补全使用信号量控制的操作,请一定不要相信某些人告诉你只要背一个什么模板&#…

Java线程入门

目录 一.线程相关概念 1.程序(program) 2.进程 3.线程 4.其他相关概念 二.线程的创建 1.继承Thread 2.Runnable接口 3.多线程机制(重要) 4.start() 三.线程终止--通知 四.线程(Thread)方法 1.常…