引言

好久好久没有更新了，来杭州的这几天真的是快把人烤熟了，好在学校的空调给力，连厕所都是中央空调，有钱真的会玩呀......

继上一篇博客

#!/usr/bin/env python

# -*- coding: utf-8 -*-

# @File : gpm_mapping.py

# @Author : zengsk in HHU

# @Time : 2019/7/21 13:40

import h5py

import numpy as np

import matplotlib as mpl

import matplotlib.pyplot as plt

from netCDF4 import Dataset as ncdataset

from mpl_toolkits.basemap import Basemap

####### 绘制GPM日降水分布图 ###########

'''

自定义色带(可替换 map.contourf()函数的cmap参数)

ps：我自己定义的实在难看，还是用自带的吧^-^

'''

def colormap():

ColorList = ['#DDDDFF','#7D7DFF','#0000C6','#000079','#CEFFCE',

'#28FF28','#007500','#FFFF93','#8C8C00','#FFB5B5',

'#FF0000','#CE0000','#750000']

return mpl.colors.LinearSegmentedColormap.from_list('cmap', ColorList, 256)

filename = r'../assets/3B-DAY.MS.MRG.3IMERG.20001026-S000000-E235959.V06.nc4'

# 这部分是读取小时分辨率数据的代码

# fHander = h5py.File(filename, mode='r')

# precip = fHander['/Grid/precipitationCal'][:]

# precip = np.transpose(precip * 0.5)

# lat = fHander['/Grid/lat'][:]

# lon = fHander['/Grid/lon'][:]

dst = ncdataset(filename)

precip = dst.variables['precipitationCal'][:]

precip = precip.reshape(3600, 1800)

precip = np.transpose(precip)

precip[np.isnan(precip)] = -999

lat = dst.variables['lat'][:]

lon = dst.variables['lon'][:]

# Plot the figure, define the geographic bounds

fig = plt.figure(3, dpi=300)

map = Basemap()

# Draw coastlines, state and country boundaries, edge of map.

map.drawcoastlines(linewidth=1.2, linestyle='solid', color='k', antialiased=3)

map.drawcountries(linewidth=0.8, linestyle='solid', color='k', antialiased=3)

# map.bluemarble(scale=3) # 加载蓝色大理石背景

# map.shadedrelief(scale=2) # 加载地形阴影

# map.etopo(scale=3) # 加载地形包括海洋

# Draw filled contours.

clevs = np.arange(0, 100, 10)

# color scale changing

# clevs = [0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0]

# Define the latitude and longitude data

x, y = np.float32(np.meshgrid(lon, lat))

# Mask the values less than 0 because there is no data to plot.

masked_array = np.ma.masked_where(precip < 0, precip)

# Plot every masked value as white

# colormap

cmap = plt.cm.get_cmap('jet')

cmap.set_bad('w', .8)

# Plot the data

cs = map.contourf(x, y, precip, clevs, cmap=cmap, latlon=True)

parallels = np.arange(-60., 61, 20.)

map.drawparallels(parallels, labels=[True, False, True, False])

meridians = np.arange(-180., 180., 60.)

map.drawmeridians(meridians, labels=[False, False, False, True])

# 设置标题和字体

plt.title('GPM IMERG Daily Rainfall( mm )')

font = {'weight': 'bold', 'size': 6}

plt.rc('font', **font)

# 加载色带

cbar = map.colorbar(cs, location='right', pad="5%")

cbar.set_label('mm/daily')

plt.savefig('../output/gpm.jpg', dpi=300) # change to your directory

plt.show()

print(" ###### 数据处理完成 ###### ")

结果展示：

image.png