python把hdf转为tif_【亲测】Python-GDAL 把MODIS HDF影像写入并导出为GeoTiff

写在前边，昨日鏖战许久未能在Windows平台上搞定基于R语言(主要是rgdal 与gdalUtils)的HDF文件的预处理，主要报错为：

Error in gdal_chooseInstallation(hasDrivers = of) :

No installations match.

然后，夜以继日地在网上找解决方案，发现只有问题，没有方案。。。内心崩溃之余，顺手换了平台与语言，以下是比较成功的完整的解决流程。

目标问题

将MODIS HDF 文件转化为R等数据分析语言可以识别并批量处理的GeoTiff格式。

问题分析

Windows + R已无法解决且网上没有比较成熟的Solution

Mac及GDAL对Python的兼容性具高

具体解决流程

1. Mac 安装 GDAL库

brew install gdal

电脑上没有Homebrew的可以看下文:

Mac 上安装Homebrew 教程

2. Python gdal模块

conda install gdal

电脑上没有Conda的可以看下文:

Mac 上安装Conda教程

3. 基于Python-GDAL将HDF文件转为GeoTiff格式

import gdal,osr

###

def array2raster(newRasterfn, rasterOrigin, xsize, ysize, array):

"""

newRasterfn: 输出tif路径

rasterOrigin: 原始栅格数据Extent

xsize: x方向像元大小

ysize: y方向像元大小

array: 计算后的栅格数据

"""

cols = array.shape[1] # 矩阵列数

rows = array.shape[0] # 矩阵行数

originX = rasterOrigin[0] # 起始像元经度

originY = rasterOrigin[1] # 起始像元纬度

driver = gdal.GetDriverByName('GTiff')

outRaster = driver.Create(newRasterfn, cols, rows, 1, gdal.GDT_Float32)

# 括号中两个0表示起始像元的行列号从(0,0)开始

outRaster.SetGeoTransform((originX, xsize, 0, originY, 0, ysize))

# 获取数据集第一个波段，是从1开始，不是从0开始

outband = outRaster.GetRasterBand(1)

outband.WriteArray(array)

outRasterSRS = osr.SpatialReference()

# 代码4326表示WGS84坐标

outRasterSRS.ImportFromEPSG(4326)

outRaster.SetProjection(outRasterSRS.ExportToWkt())

outband.FlushCache()

###

from os import listdir

import re

inputf = '/Users/jerseyshen/Documents/Global_NDVI/'

output = '/Users/jerseyshen/Documents/Global_NDVI_new/'

files = listdir(inputf)

pattern = '.hdf$'

timer = 0

for i in files:

temp_input = inputf+i

temp_out_name = re.split(pattern,i)[0]+'.tif'

temp_output = output+temp_out_name

df = ds = gdal.Open(temp_input)

subdatasets = ds.GetSubDatasets()

ndvi_ds = gdal.Open(subdatasets[0][0]).ReadAsArray()

xsize = 0.05

ysize = 0.05

#对于全球尺度rasterOrigin 为[0(xmin),-90(ymin)]

array2raster(temp_output,[0,-90],xsize,ysize,ndvi_ds)

timer = timer + 1

print(timer)

4. GeoTiff文件在R语言中的呈现

image.png