概述（前言）

本文作为上一篇博客《MATLAB处理GPM(IMERG/GSMaP)卫星降水数据》的后续，作者选择另一种语言python来实现处理降水领域高分辨率降水数据GPM-IMERG和GPM-GSMaP，两种数据的详细介绍、下载方式以及matlab读取参见上一篇文章《MATLAB处理GPM(IMERG/GSMaP)卫星降水数据》：

image.png

正文

Python 处理GPM-IMERG数据

1. 数据下载(python ftplib)
参见文章《python 实现FTP文件批量下载》
2.数据读取(python)

#!/usr/bin/env python
# -*- coding: utf-8 -*-
# @Time    : 2019/2/28 17:28
# @Author  : zengsk in HoHai
# @File    : Calc_AVG_IMERG.py

'''
Calculating the average precipitation over a given region
'''

# Import packages
import os
import h5py
import numpy as np
import pandas as pd
import warnings

warnings.simplefilter("ignore")

# header
header = '''ncols         3600
nrows         1800
xllcorner     -180
yllcorner     -90
cellsize      0.1
NODATA_value  -999.0'''

# dataset directory
sPath = r'D:\降水数据\satellite_pre\IMERG\late\201508'
# mask region
mask = pd.read_table(r'.\assets\china_mask_gpm_01.txt', sep='\s+', header=None, skiprows=6).values

rain = np.empty([1800, 3600], dtype=float)
tag = np.empty(rain.shape, dtype=float)     # Record the effective precip times per grid

Files = os.listdir(sPath)
for file in Files:
   fpath = os.path.join(sPath, file)
   if os.path.splitext(fpath)[1] == '.RT-H5':
       print(os.path.basename(fpath))
       f = h5py.File(fpath)      # Open the H5 file, return the File class
       Uncal = f['/Grid/precipitationUncal'].value

       Uncal = np.transpose(Uncal)   # size: 1800 x 3600
       Uncal = np.flipud(Uncal)      # convert to  90N ~ 90S
       tag[Uncal>=0] = tag[Uncal>=0] + 1
       Uncal[Uncal<0] = 0
       Uncal /= 2
       rain += Uncal

# output
rainfall_avg = rain / tag * 2 * 24    # unit conversion
rainfall_avg[np.isnan(rainfall_avg)] = -999.00 # NODATA_value
rainfall_avg[mask<0] = -999.00

ofname = r'.\output\IMERG_AVG_201508mon.txt'
np.savetxt(ofname, rainfall_avg, fmt='%7.2f ', comments='', header=header)

print("\n***************** 数据处理完成!!! Nice Job!!! ****************")

---

Python 处理GPM-GSMaP数据

1. 数据下载(python ftplib)
参见文章《python 实现FTP文件批量下载》
2.数据读取(python)

#!/usr/bin/env python
# -*- coding: utf-8 -*-
# @Time    : 2019/4/13 18:34
# @Author  : zengsk in HoHai
# @File    : Calc_Avg_GSMaP.py

'''
Calculating the average precipitation over a given region
 GsMaP info:
    60N-->60S, 0-->360E;
    0.1*0.1degree;
    1 hourly or daily;
    4-byte-float;
'''

import os
import struct
import numpy as np
import pandas as pd
import warnings

warnings.simplefilter("ignore")

# header
header = '''ncols         3600
nrows         1200
xllcorner     0
yllcorner     -60
cellsize      0.1
NODATA_value  -999.0'''

# dataset directory
sPath = r'D:\data\GSMaP'
# mask region file
mask = pd.read_table(r'.\assets\china_mask_gpm_01.txt', sep='\s+', header=None, skiprows=6).values

rain = np.empty([1200, 3600], dtype=float)
tag = np.empty(rain.shape, dtype=float)     # Record the effective precip times per grid

Files = os.listdir(sPath)
for file in Files:
    filepath = os.path.join(sPath, file)
    if '.dat' == os.path.splitext(filepath)[1]:
        print(filepath)
        with open(filepath, 'rb') as fid:
            data = fid.read()
            data = struct.unpack('f'*3600*1200, data)
        data = np.array(data)
        data.resize(1200,3600)
        tag[data >= 0] = tag[data >= 0] + 1
        data[data < 0] = 0.0
        data *= 24           # unit conversion
        rain += data
# output
rainfall_avg = rain / tag
rainfall_avg[np.isnan(rainfall_avg)] = -999.00 # NODATA_value
rainfall_avg[mask<0] = -999.00

ofname = r'.\output\201607avg_gsmap.txt'
np.savetxt(ofname, rainfall_avg, fmt='%7.2f ', comments='', header=header)

print("\n+---------------- 数据处理完成!!! Nice Job!!! -----------------+\n")

---

最后来一张IMERG 2015年08月全球平均降水量分布图(由上述代码读取)

IMERG 2015年08月全球平均降水量分布

欢迎留言讨论！转载请注明出处！！！