【Python实例】Python读取并绘制nc数据
- 准备:安装netCDF库等
- 读取nc数据相关信息
- 绘制图形
- 利用basemap绘图
- 参考
准备:安装netCDF库等
以【1960-2020年中国1km分辨率月降水数据集】中2020年降水为例。
先在Panopoly中查看数据属性,如下:
右键绘制图形,如下:
读取nc数据相关信息
Python代码如下:
import netCDF4
from netCDF4 import Dataset
import matplotlib.pyplot as plt
import numpy as np
# 打开下载的.nc文件
ncfile = Dataset("D:/0 DataBase/CN_Prec_MonthlyMean_1km_2020.nc/CN_Prec_MonthlyMean_1km_2020.nc")
print(ncfile.variables.keys())
# 读取变量
lat = ncfile.variables['latitude'][:]
long = ncfile.variables['longitude'][:]
[X, Y] = np.meshgrid(lat, long) # 合成矩阵,将向量latitude和longitude定义的区域转换成矩阵X和Y
pre = ncfile.variables['pre'][:]
time = ncfile.variables['time'][:]
# 绘制图形
plt.contourf(Y, X,pre[0])
plt.colorbar(label="MonthlyMean Precipitation", orientation="horizontal")
plt.show()
# 关闭nc文件
ncfile.close()
输出属性如下:
绘制图形
利用basemap绘图
Python代码如下:
#利用basemap绘图
import netCDF4
from mpl_toolkits import basemap
from netCDF4 import Dataset
from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np
# 打开下载的.nc文件
ncfile = Dataset("D:/0 DataBase/CN_Prec_MonthlyMean_1km_2020.nc/CN_Prec_MonthlyMean_1km_2020.nc")
# 读取变量
lat = ncfile.variables['latitude'][:]
long = ncfile.variables['longitude'][:]
[X, Y] = np.meshgrid(lat, long) # 合成矩阵,将向量latitude和longitude定义的区域转换成矩阵X和Y
pre = ncfile.variables['pre'][:]
time = ncfile.variables['time'][:]
# 创建一个图形和一个轴
fig = plt.figure(figsize=(10, 8))
ax = fig.add_subplot(111)
# 设定标题
titleName = "MonthlyMean Precipitation"
# 初始化Basemap对象,设置地图参数
# 设置投影方式:cyl为圆柱投影、还可设置merc为mercator投影,lcc兰伯特投影等。
# 参数 "resolution" 用于控制地图面积边缘的精细程度,有 'crude', 'low' and'intermediate' resolution datasets are installed by default.
# 如果需要'high' or 'full' resolution dataset,需要`conda install -c conda-forge basemap-data-hires`.
m = Basemap(lat_1=30, lat_2=62, lon_0=103, projection='lcc', resolution='i',
llcrnrlon=80, llcrnrlat=10, urcrnrlon=145, urcrnrlat=55)
# 设置经纬度
# 绘制纬度坐标,纬度每隔10度画一条线,且左侧和右侧标注纬度。
m.drawparallels(np.arange(10., 55., 10.), labels=[1, 1, 0, 0], fontsize=10)
# 绘制经度坐标,经度每隔10度画一条线,且上侧和下侧标注经度。
m.drawmeridians(np.arange(70., 160., 10.), labels=[0, 0, 1, 1], fontsize=10)
# 形状文件的路径
#shp_path = "D:/0 DataBase/0 GBA Database/Boundary_Projected"
# 添加中国地图图层
#m.readshapefile(shp_path,'Boundary_Projected', drawbounds=True)
# 添加大陆地图图层
m.fillcontinents(color="coral", lake_color="blue")
# 设定海洋的颜色
m.drawmapboundary(fill_color="aqua")
# 在指定位置绘制坐标点
lon_me = long.mean()
lat_me = lat.mean()
lon_do, lat_do = np.meshgrid(lon_me, lat_me) # 生成网格点坐标矩阵
# 转换经纬度到Basemap坐标系
xi, yi = m(lon_do, lat_do)
# print(xi,yi)
m.plot(xi, yi, "ro") # 绘制红色点
# 在指定位置标注文字
plt.text(xi, yi, "A1(%5.1fW, %3.1fN)" % (lon_me, lat_me))
# 关闭.nc文件
ncfile.close()
# 显示图形
plt.show()
图形绘制如下: