| 从xarray走向netCDF处理(四)：合并与计算 | 您所在的位置：网站首页 › 经纬度加减乘除怎么计算 › 从xarray走向netCDF处理(四)：合并与计算 |
从xarray走向netCDF处理(四)：合并与计算
|
前面有关xarray已经讲了3期了，介绍了数据索引，数据结构还有插值和掩膜。今天这是最后一期介绍用xarray处理nc数据了，打算聊一下如何做数据合并与计算。 数据合并数据合并主要是两种形式 维度的拼接：如将日数据合成为年数据，就属于在时间维度上的合并。变量的合并：如将多个物理量合到同一个Dataset中。xarray围绕着这两种合并方式介绍了concatenate, merge, combine, update四种方法。 我在这里就挑最常用的跟大家聊聊。 维度拼接使用 concat() 方法可以实现维度的拼接。 下面是演示数据，来源于2018年和2019年前三个月的ERA-Interim月平均数据。 1>>> ds2018 2 3Dimensions: (latitude: 241, longitude: 480, time: 12) 4Coordinates: 5 * longitude (longitude) float32 0.0 0.75 1.5 2.25 ... 357.75 358.5 359.25 6 * latitude (latitude) float32 90.0 89.25 88.5 87.75 ... -88.5 -89.25 -90.0 7 * time (time) datetime64[ns] 2018-01-01 2018-02-01 ... 2018-12-01 8Data variables: 9 u10 (time, latitude, longitude) float32 ... 10 v10 (time, latitude, longitude) float32 ... 11 t2m (time, latitude, longitude) float32 ... 12Attributes: 13 Conventions: CF-1.6 14 15>>> ds2019 16 17Dimensions: (latitude: 241, longitude: 480, time: 3) 18Coordinates: 19 * longitude (longitude) float32 0.0 0.75 1.5 2.25 ... 357.75 358.5 359.25 20 * latitude (latitude) float32 90.0 89.25 88.5 87.75 ... -88.5 -89.25 -90.0 21 * time (time) datetime64[ns] 2019-01-01 2019-02-01 2019-03-01 22Data variables: 23 u10 (time, latitude, longitude) float32 ... 24 v10 (time, latitude, longitude) float32 ... 25 t2m (time, latitude, longitude) float32 ... 26Attributes: 27 Conventions: CF-1.6ds2018时间维度为12，ds2019时间维度为3，下面使用 concat() 合并后时间维度为15 1>>> xr.concat([ds2018, ds2019], dim='time') 2 3Dimensions: (latitude: 241, longitude: 480, time: 15) 4Coordinates: 5 * longitude (longitude) float32 0.0 0.75 1.5 2.25 ... 357.75 358.5 359.25 6 * latitude (latitude) float32 90.0 89.25 88.5 87.75 ... -88.5 -89.25 -90.0 7 * time (time) datetime64[ns] 2018-01-01 2018-02-01 ... 2019-03-01 8Data variables: 9 u10 (time, latitude, longitude) float32 -0.9599868 ... 4.5229325 10 v10 (time, latitude, longitude) float32 3.1737509 ... -2.289166 11 t2m (time, latitude, longitude) float32 248.46857 ... 225.19632 12Attributes: 13 Conventions: CF-1.6变量合并使用 merge() 方法，可以将ds2018中的u10和ds2019中的t2m合并到一起，而且在时间维上缺失会自动设置为nan。 1>>> xr.merge([ds2018.u10, ds2019.t2m]) 2 3Dimensions: (latitude: 241, longitude: 480, time: 15) 4Coordinates: 5 * time (time) datetime64[ns] 2018-01-01 2018-02-01 ... 2019-03-01 6 * longitude (longitude) float32 0.0 0.75 1.5 2.25 ... 357.75 358.5 359.25 7 * latitude (latitude) float32 90.0 89.25 88.5 87.75 ... -88.5 -89.25 -90.0 8Data variables: 9 u10 (time, latitude, longitude) float32 -0.9599868 -0.9599868 ... nan 10 t2m (time, latitude, longitude) float32 nan nan ... 225.19632数据计算最基本的计算就是进行加减乘除，任意一个DataArray或者Dataset都可以直接进行四则运算。 除此以外，xarray还可以帮你快速地求出平均值，方差，最小值，最大值等。你可以指定具体对那个维度进行计算，如果不指定维度默认会对所有维度进行计算。 比如要对经、纬两个维度进行平均，最后的结果只有时间维的12个值。 1>>> ds2018.mean(dim=['latitude', 'longitude']) 2 3Dimensions: (time: 12) 4Coordinates: 5 * time (time) datetime64[ns] 2018-01-01 2018-02-01 ... 2018-12-01 6Data variables: 7 u10 (time) float32 -0.120867714 0.13738841 ... -0.016953295 0.08418254 8 v10 (time) float32 -0.21417202 -0.12471106 ... -0.082666814 9 t2m (time) float32 277.57446 277.32916 276.72095 ... 278.1613 278.01758如果对时间维进行求方差，则结果会保留空间场上的信息。 1>>> ds2018.std(dim='time') 2 3Dimensions: (latitude: 241, longitude: 480) 4Coordinates: 5 * longitude (longitude) float32 0.0 0.75 1.5 2.25 ... 357.75 358.5 359.25 6 * latitude (latitude) float32 90.0 89.25 88.5 87.75 ... -88.5 -89.25 -90.0 7Data variables: 8 u10 (latitude, longitude) float32 1.9192954 1.9192954 ... 1.2133 9 v10 (latitude, longitude) float32 1.3066719 1.3066719 ... 1.577495 10 t2m (latitude, longitude) float32 9.5681305 9.5681305 ... 11.313364而且xarray在时间维上的计算还有很多贴心的用法，比如月数据转年数据，月数据转季节数据。 1>>> ds2018.groupby('time.season').min(dim='time') 2 3Dimensions: (latitude: 241, longitude: 480, season: 4) 4Coordinates: 5 * longitude (longitude) float32 0.0 0.75 1.5 2.25 ... 357.75 358.5 359.25 6 * latitude (latitude) float32 90.0 89.25 88.5 87.75 ... -88.5 -89.25 -90.0 7 * season (season) object 'DJF' 'JJA' 'MAM' 'SON' 8Data variables: 9 u10 (season, latitude, longitude) float32 -3.2505732 ... 3.4194348 10 v10 (season, latitude, longitude) float32 0.9903783 ... -3.318876 11 t2m (season, latitude, longitude) float32 248.46857 ... 214.61685一个栗子 1import numpy as np 2import xarray as xr 3import cartopy.crs as ccrs 4import cartopy.feature as cfeat 5from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER 6import matplotlib.pyplot as plt 7 8def create_map(ax): 9 # 设置地图属性 10 ax.add_feature(cfeat.BORDERS.with_scale('50m'), linewidth=0.8) # 加载分辨率为50的国界线 11 ax.add_feature(cfeat.COASTLINE.with_scale('50m'), linewidth=0.6) # 加载分辨率为50的海岸线 12 ax.add_feature(cfeat.RIVERS.with_scale('50m')) # 加载分辨率为50的河流 13 ax.add_feature(cfeat.LAKES.with_scale('50m')) # 加载分辨率为50的湖泊 14 # 设置网格点属性 15 gl = ax.gridlines(crs=ccrs.PlateCarree(), draw_labels=True, 16 linewidth=1.2, color='k', alpha=0.5, linestyle='--') 17 gl.xlabels_top = False # 关闭顶端的经纬度标签 18 gl.ylabels_right = False # 关闭右侧的经纬度标签 19 gl.xformatter = LONGITUDE_FORMATTER # x轴设为经度的格式 20 gl.yformatter = LATITUDE_FORMATTER # y轴设为纬度的格式 21 return ax 22 23if __name__ == '__main__': 24 # 创建画图空间 25 proj = ccrs.PlateCarree() 26 fig, axes = plt.subplots(nrows=4, ncols=1, figsize=(7, 14), subplot_kw={'projection': proj}) 27 ds = xr.open_dataset('EC-Interim_monthly_2018.nc') 28 temp = (ds['t2m'] - 273.15).groupby('time.season').mean('time') 29 # --画图 30 cbar_kwargs = { 31 'label': '2m temperature (℃)', 32 'ticks': np.arange(-30, 30+5, 5) 33 } 34 levels = np.arange(-30, 30+1, 1) 35 for i, season in enumerate(('DJF', 'MAM', 'JJA', 'SON')): 36 ax = create_map(axes[i]) 37 temp.sel(season=season).plot.contourf(ax=ax, levels=levels, cmap='Spectral_r', extend='both', 38 cbar_kwargs=cbar_kwargs, transform=ccrs.PlateCarree()) 39 fig.show()
|
【本文地址】
公司简介
联系我们