3. 中尺度气象数据应用:ERA5、MERRA-2数据下载与预处理,网格点插值方法

做风资源评估这些年,我越来越觉得,中尺度气象数据就像是我们手里的“望远镜”。测风塔只能看到几十米高的风,而中尺度数据能让我们看到几十公里甚至几百公里范围内的风场演变。说白了,没有这个“望远镜”,微观选址就像在黑夜里摸路。

今天咱们就聊聊ERA5和MERRA-2这两个主流数据集。我个人的习惯是,项目前期先用它们做区域风资源普查,等确定几个候选机位点后,再用测风塔数据做精细校正。嗯,这里要注意,中尺度数据不能直接替代测风塔,但能帮我们少走很多弯路。

3.1 ERA5与MERRA-2:两个“老伙计”的对比

先说说这两个数据集的特点。ERA5是欧洲中期天气预报中心(ECMWF)的产品,MERRA-2则是美国NASA的杰作。我在项目中两个都用过,各有千秋。

特性 ERA5 MERRA-2
空间分辨率 0.25° × 0.25°(约31km) 0.5° × 0.625°(约50km)
时间分辨率 逐小时 逐小时
时间跨度 1940年至今 1980年至今
垂直层数 137层(地表到0.01hPa) 72层(地表到0.01hPa)
数据更新 滞后约3个月 滞后约2个月
风速变量 10m u/v分量,100m u/v分量 10m u/v分量,50m u/v分量

你可能会问,为什么MERRA-2有50m风速而ERA5有100m?这跟两个机构的模型设计有关。我个人更偏爱ERA5,因为它的空间分辨率更高,而且直接提供了100m高度的风速——这对现代大型风电机组来说太重要了。不过MERRA-2也有优势,它的同化系统对气溶胶和地表过程的处理更精细,在复杂地形区域表现反而更好。

关键点:ERA5适合平坦地形和海上风电项目,MERRA-2在复杂山地和干旱地区有独特优势。我建议两个都下载,交叉验证。

3.2 数据下载:别让第一步卡住你

下载数据听起来简单,但坑不少。我曾经因为没处理好API密钥,折腾了一整天。来,直接上干货。

3.2.1 ERA5下载(CDS API)

ERA5通过CDS(Climate Data Store)提供数据。你需要先注册账号,获取API密钥。嗯,这一步很多人会忽略——记得把密钥保存到 .cdsapirc 文件中。

# 安装CDS API
pip install cdsapi

# Python下载脚本示例
import cdsapi

c = cdsapi.Client()

c.retrieve(
    'reanalysis-era5-single-levels',
    {
        'product_type': 'reanalysis',
        'format': 'netcdf',
        'variable': [
            '10m_u_component_of_wind', '10m_v_component_of_wind',
            '100m_u_component_of_wind', '100m_v_component_of_wind',
        ],
        'year': '2023',
        'month': ['01', '02', '03', '04', '05', '06',
                  '07', '08', '09', '10', '11', '12'],
        'day': [f'{d:02d}' for d in range(1, 32)],
        'time': [f'{h:02d}:00' for h in range(24)],
        'area': [40, 110, 20, 125],  # 北纬40°到20°,东经110°到125°
    },
    'era5_wind_2023.nc')

小技巧:下载时尽量选择NetCDF格式,它比GRIB更容易处理。另外,一次下载一整年的数据,别一个月一个月地下载,否则会被CDS限流。

3.2.2 MERRA-2下载

MERRA-2通过NASA GES DISC提供。下载方式稍微不同,需要用 earthaccess 库。

# 安装earthaccess
pip install earthaccess

# Python下载脚本
import earthaccess
import xarray as xr

# 登录NASA Earthdata
earthaccess.login(strategy="interactive")

# 搜索数据
results = earthaccess.search_data(
    short_name='M2T1NXSLV',  # 单层变量
    bounding_box=(110, 20, 125, 40),
    temporal=('2023-01-01', '2023-12-31'),
    count=10
)

# 下载
files = earthaccess.download(results, './merra2_data/')

注意:NASA Earthdata需要单独注册账号。我遇到过好几次下载到一半断连的情况,建议用 earthaccess 的断点续传功能。

3.3 数据预处理:从原始数据到可用格式

下载下来的数据不能直接用。你想想看,原始NetCDF文件里包含了全球数据,而我们只需要项目区域的那一小块。另外,风速是以u/v分量存储的,需要合成标量风速和风向。

3.3.1 区域裁剪与变量提取

import xarray as xr
import numpy as np

# 读取ERA5数据
ds = xr.open_dataset('era5_wind_2023.nc')

# 裁剪到项目区域(例如:北纬35°-38°,东经115°-118°)
ds_region = ds.sel(
    latitude=slice(38, 35),  # 注意纬度从大到小
    longitude=slice(115, 118)
)

# 提取100m风速分量
u100 = ds_region['u100']
v100 = ds_region['v100']

# 合成风速和风向
ws100 = np.sqrt(u100**2 + v100**2)
wd100 = np.arctan2(v100, u100) * 180 / np.pi
wd100 = (270 - wd100) % 360  # 转换为气象风向

# 保存为新的NetCDF
ds_out = xr.Dataset({
    'ws100': ws100,
    'wd100': wd100
})
ds_out.to_netcdf('era5_region_ws_wd.nc')

避坑指南:我曾经在风向转换上栽过跟头。ERA5的u/v分量是气象坐标系(u为东向正,v为北向正),而我们需要的是“风从哪里来”的气象风向。公式是:wd = (270 - arctan2(v, u) * 180/π) % 360。别搞反了!

3.3.2 时间序列质量控制

中尺度数据偶尔会有异常值。我一般会做三步检查:

  1. 范围检查:风速0-50m/s,风向0-360°
  2. 时序连续性:相邻小时风速变化不超过20m/s
  3. 气候学检查:与长期平均值偏差不超过3个标准差
# 简单的质量控制
def quality_control(ws, wd):
    # 标记异常值
    mask = (ws < 0) | (ws > 50) | (wd < 0) | (wd > 360)
    
    # 标记时序跳变
    ws_diff = np.abs(np.diff(ws, axis=0))
    mask[1:] |= (ws_diff > 20)
    
    # 用前后平均值插补
    ws_clean = ws.copy()
    ws_clean[mask] = np.nan
    ws_clean = ws_clean.interpolate_na(dim='time', method='linear')
    
    return ws_clean

3.4 网格点插值方法:把“粗格子”变“细”

中尺度数据的网格分辨率是几十公里,而我们的机位点可能就在两个网格点之间。这时候就需要插值。说白了,就是根据周围几个网格点的值,估算出目标点的值。

我常用的插值方法有三种,各有适用场景:

方法 原理 适用场景 精度
双线性插值 用周围4个网格点线性加权 平坦地形,快速计算 中等
反距离加权(IDW) 距离越近权重越大 复杂地形,局部效应明显 较高
克里金插值 考虑空间自相关性 需要统计最优估计 最高

3.4.1 双线性插值(最常用)

from scipy.interpolate import RegularGridInterpolator

# 假设目标点经纬度
target_lon, target_lat = 116.5, 36.8

# 获取网格坐标
lons = ds_region['longitude'].values
lats = ds_region['latitude'].values

# 创建插值函数
interp_func = RegularGridInterpolator(
    (lats, lons), 
    ws100.values,
    method='linear'
)

# 插值到目标点
ws_target = interp_func((target_lat, target_lon))

经验之谈:双线性插值在平坦地形上误差很小,但在山地,我建议用IDW。因为山地风场受地形影响大,距离近的网格点更有参考价值。

3.4.2 反距离加权插值(IDW)

from scipy.spatial import cKDTree

def idw_interpolation(lons, lats, values, target_lon, target_lat, power=2):
    # 构建KD树
    points = np.column_stack([lons.ravel(), lats.ravel()])
    values_flat = values.ravel()
    
    tree = cKDTree(points)
    
    # 查找最近的4个点
    distances, indices = tree.query([target_lon, target_lat], k=4)
    
    # 计算权重
    weights = 1.0 / (distances ** power)
    weights /= weights.sum()
    
    # 加权平均
    return np.sum(values_flat[indices] * weights)

# 使用
ws_target = idw_interpolation(
    lons, lats, ws100.values,
    target_lon, target_lat
)

3.5 实战案例:从数据到风资源图谱

我记得去年在内蒙古做一个项目,项目区域有200平方公里,但只有一座测风塔。我用了ERA5数据,做了以下几步:

  1. 下载项目区域及周边2°×2°范围的ERA5逐小时数据
  2. 用双线性插值生成1km×1km的网格风场
  3. 用测风塔数据做偏差校正(线性回归)
  4. 生成年平均风速分布图和风功率密度图

结果发现,测风塔位置的风速是7.2m/s,但项目区域东北角有个小山脊,插值结果显示风速能达到8.1m/s。后来实地踏勘,确实在那找到了一个很好的机位点。这就是中尺度数据的价值——它能帮你发现测风塔“看不到”的风资源。

核心要点:中尺度数据是“面”上的工具,测风塔是“点”上的工具。两者结合,才能做出靠谱的微观选址。

3.6 知识体系总览

下面这张图是我自己整理的,把这一章的核心逻辑串起来了。你一看就明白:从数据源到预处理,再到插值方法,最后到实际应用,每一步都有坑,也都有技巧。

中尺度气象数据应用流程 数据源 数据源 ERA5 (0.25° × 0.25°) MERRA-2 (0.5° × 0.625°) 数据预处理 区域裁剪 u/v → 风速/风向 质量控制 网格点插值方法 双线性插值 反距离加权 (IDW) 克里金插值 风资源图谱生成 · 微观选址优化

嗯,这一章的内容就这些。数据下载、预处理、插值,每一步都有门道。我建议你拿到数据后,先花点时间做质量控制,别急着往下跑。数据质量决定了最终结果的可靠性——这个道理,我在无数个项目里验证过了。

专注资料整理