ArcPy常用函数示例01-栅格数据投影转换

我不是证明我有多了不起，我是要证明我失去的东西我一定要亲手拿回来（“英雄本色”-周润发）。

---

前言

  1. 概述

• 基于arcpy.ProjectRaster_management函数对栅格数据进行投影转换
• 如需程序测试数据，请留下邮箱

  2. 版本

   2.1 山东青岛，2021年5月20日，Version 1

  3. 参考资料

[1]. ArcPy函数帮助手册

---

一、投影转换的目的

1. 地理数据必须有地理坐标参考或者投影坐标参考，不同坐标参考下，对于同一栅格数据，其栅格形状、数据的行数和列数不一样，当我们需要依据栅格行列号进行不同影像间的对应分析时，就必须确保不同影像之间坐标系统的一致性，才能做到不同影像间像元的一一对应。
2. 不同投影坐标系，侧重不同，比如当我们需要进行面积分析时，就需要 “等积投影”。

---

二、arcpy.ProjectRaster_management()函数

1. 函数位置
   借助ArcPy自带函数实现栅格数据投影转换，下图展示了该函数的位置

2. 函数用法
   2.1 函数的输入参数
      该函数一共8个输入参数，其中三个必须输入参数，分别是
   （1）in_raster，输入栅格数据，即需要转换投影的栅格数据
   （2）out_raster，输出栅格数据，即转换投影之后的栅格数据
   （3）out_coor_system，输出坐标系统，即想要转换成的坐标系统

   2.2 说明
       需要注意resampling_type即重采样方法的选择。 NEARST是默认的重采样方法，该方法和“Majourity”方法适合于“categorical data”即分类数据。对于连续数据建议使用BILINEAR方法或者CUBIC方法。

---

三、应用示例

1. 测试代码
   测试代码的功能是将某一文件夹下的所有.tif数据进行投影转换，由地理坐标系转换为等积投影，测试代码如下：

# -*- coding: cp936 -*-

'''
1. 程序目的
  (1) 将GCS-WGS1984地理坐标系转换为Asia_North_Albers_Equal_Area_Conic.prj投影坐标系

2. 山东青岛  2021年5月20日

3. 说明
   (1) 程序调用ArcPy的函数，所以是在ArcMap自带的Python27环境下编写的程序
'''

# 0. 相关包的导入
import os
import glob

import arcpy;from arcpy import env
from arcpy.sa import *
import shutil
import datetime
arcpy.CheckOutExtension("Spatial")

# 1. 路径处理和基本变量定义

  # (1) 输入路径
inpath = '存储.tif数据的文件夹完整路径\\'

  # (2) 输出路径
outpath = inpath + 'Project'
if os.path.exists(outpath):
    shutil.rmtree(outpath)
    os.mkdir(outpath)
else:
    os.mkdir(outpath)

  # (3) 目标投影坐标系   
EqualAreaPro = rootdir + 'Asia_North_Albers_Equal_Area_Conic.prj'
  # (4) 重新投影之后的像元大小
cellsize = 50000 # 50km

# 2. 调用arcpy.ProjectRaster_management()函数，实现批量投影转换
  # (1) 获取某一路径下的所有.tif文件
in_raster_all = glob.glob(os.path.join('%s%s')%(inpath,'*.tif'))
  # (2) 循环批量投影转换
for ii in range(len(in_raster_all)):
    in_raster = in_raster_all[ii] # 输栅格数据
    dir_path,inraster_name = os.path.split(in_raster) # 分离路径和输入栅格名称
    outraster_name = inraster_name.replace('.tif','_Albers-EA.tif')
    out_raster = outpath + '\\' + outraster_name

    # 批量投影转换
    arcpy.ProjectRaster_management(in_raster,out_raster,EqualAreaPro,"NEAREST",cellsize,\
                            "#", "#", "#")
    
    print outraster_name

print '*******************************************************'
print 'Finished'

---