Python讀取GSMap數據的問題

glob:是一個文件操作相關模塊，支持通配符操作，可用於查找符合一定條件的文件，這裡用來批量查找".dat"文件
numpy:讀取“.dat”文件為array，固定尺寸為[1200, 3600]，填補異常值和缺失值為-999.0
arcpy:將numpy的array輸出為“.tiff”格式，添加空間參考GCS_WGS_1984，最終保存tif文件

glob.glob:用於選擇“.dat”格式文件
np.fromfile:讀取二進制文件到array
arcpy.NumPyArrayToRaster:將array轉化為柵格數據
arcpy.SpatialReference(4326):引用GCS_WGS_1984的空間參考代號
arcpy.DefineProjection_management:建立空間參考

#!/usr/bin/env python
# -*- coding:utf-8 -*-
# Name : gsmap2TIFF.py
# Author : zengsk in NanJing
# Created: 2019/8/24 23:26

'''
說明：1.該腳本是讀取gsmap小時尺度降水數據，輸出為tiff
 2.運行環境需要安裝python2 需要arcpy模塊
 3.使用arcgis自帶的python環境(有arcpy模塊)
 4.運行結果可以直接用arcgis打開
'''

# 導入module
import os
import glob
import numpy as np
import arcpy
import warnings

warnings.simplefilter("ignore") # 忽略警告

# 原始降水數據文件夾,可根據自己本地情況修改
sPath = r'./data'
oDir = r"./tiff/"

for fileName in glob.glob(sPath + '\*.dat'):
 print("Processing... {0}".format(fileName))
 ds = np.fromfile(fileName, dtype=np.float32)
 ds = np.resize(ds, [1200, 3600])
 ds[ds < 0] = -999.00
 ds[np.isnan(ds)] = -999.00 # NODATA_value

 # 輸出為TIFF(註意：要用到arcpy模塊)
 if not os.path.exists(oDir):
  os.makedirs(oDir)
 TiffName = oDir + os.sep + os.path.basename(fileName)[0:-4] + '.tif' # 輸出文件名(可根據實際情況改)
 # arcpy.NumPyArrayToRaster()不清楚輸入參數可以查看arcpy的官方文檔
 # 矩陣轉為柵格
 raster = arcpy.NumPyArrayToRaster(ds, arcpy.Point(0, -60.0),
          x_cell_size=0.1, y_cell_size=0.1, value_to_nodata=-999.00)
 # 添加地理坐標系 GCS_WGS_1984
 spatialRef = arcpy.SpatialReference(4326)
 arcpy.DefineProjection_management(raster, spatialRef)
 raster.save(TiffName)

print("\n++++++ Data Processing Successfully Completed ! ++++++")