【python/MIKE21】使用mikeio库制作风场输入dfs2文件--以era5风场为例。

该文参考博主@Anywayyyyy.的博文《如何使用python处理nc数据制作Mike风场文件–以ERA5数据为例》
链接: link
代码结构主体参考上文，但由于@Anywayyyyy.博主使用的mikeio库是1.0以下版本，1.0以上版本进行了大量改动导致无法直接使用，而安装老版本的mikeio库也会出现各种问题，因此笔者查找了源码及注释，给出了新版本mikeio库制作风场dfs2文件的方法，代码如下。

"""
把era5下载的风场改写为mike的输入文件dfs2格式
"""
import mikeio
from mikeio import Dfs2
from mikeio.eum import ItemInfo, EUMUnit
from mikeio.eum import EUMType
import numpy as np
import netCDF4 as nc
import datetime as dt

file = nc.Dataset('D:/MIKE Zero Projects/waigaoqiao/era5_wind_and_pressure.nc', 'r')  # 文件路径及名称
time = np.array(file.variables['time'][:])  # 读取时间信息
lon = np.array(file.variables['longitude'][:])  # 读取经度信息
lat = np.flipud(np.array(file.variables['latitude'][:]))  # 读取纬度信息
u = np.flip(np.array(file.variables['u10'][:]), axis=1)  # 读取东西风速u
v = np.flip(np.array(file.variables['v10'][:]), axis=1)  # 读取南北风速v
p = np.flip(np.array(file.variables['sp'][:]), axis=1)  # 读取表面压强

由于era5风场数据的纬度是从大到小排列的，与平时的习惯不同，为了防止后面出现各种问题，在导入数据后将其以及其他变量进行翻转以符合正常习惯。

u = u * 1.1  # 将风速放大1.1倍
v = v * 1.1

time_dt = []  # 指定一个新的列表
tstart = dt.datetime(1900, 1, 1, 0)  # nc文件的起始时刻
for i in time:
    time_dt.append(tstart+dt.timedelta(hours=int(i)))  # 将nc文件中的时刻转换为正常的时间格式

dfsfilename = 'D:/MIKE Zero Projects/waigaoqiao/2021uvp_for_fm.dfs2'  # 给定文件名
x0 = lon[0]
y0 = lat[0]
dx = 0.25  # 设置网格间距，ERA5默认网格间距均为0.25度
dy = 0.25
time0 = time_dt[0]   # 设置其实时间，用到了上一步的时间处理
delta_t = 3600   # 设置时间间隔，ERA5数据时间间隔为1h，3600s

# 创建变量， EUMType.Wind_Velocity说明变量类型为风速, EUMUnit.meter_per_sec说明变量单位为m/s
items = [ItemInfo("u", EUMType.Wind_Velocity, EUMUnit.meter_per_sec),
         ItemInfo("v", EUMType.Wind_Velocity, EUMUnit.meter_per_sec),
         ItemInfo("p", EUMType.Pressure, EUMUnit.pascal)]

geometry = mikeio.Grid2D(x0=x0, nx=len(lon), dx=dx, y0=y0, ny=len(lat), dy=dy, projection="LONG/LAT")

dfs = Dfs2()  # 实例化一个dfs变量

d = [np.array(u), np.array(v), np.array(p)]  # 将要输出至dfs2的数组打包
aa = mikeio.Dataset(data=d, time=time_dt, items=items, geometry=geometry)

# dfs.write(filename=dfsfilename, data=d, start_time=time0, dt=delta_t, items=items, coordinate=coordinate, dx=dx, dy=dy)  # 输出变量
dfs.write(data=aa, filename=dfsfilename, dt=delta_t, dx=dx, dy=dy, title='wind_uvp')  # 输出变量

写的比较粗糙，有各种问题欢迎大家评论区留言，笔者会根据问题再做修改。