跑过WRF-GHG、又碰过VPRM通量的人都有同感气象场上手容易天文台的GFS/FNL打包好就能喂给WPS但温室气体浓度边界这一关稍不留神就能让整个模拟域的CO₂、CH₄背景场跑偏。前几篇文章把WRF-GHG-Prepy工具从目录结构讲到气候态背景场处理这篇继续往下走专门拆CAMS-Inversion数据这条线。文章会先把CAMS-Inversion在边界条件里的定位说清楚再带着把Prepy里对应的处理函数一层层剥开最后给一套可直接照抄的实操流程和排错清单。适合谁看正在跑WRF-GHG或者说跑WRF-VPRM耦合模拟的同行。如果你已经配好了前置工具手头有WPS生成的met_em文件这篇能帮你把全球反演浓度场正确插进WRF的侧边界和初始场里。如果还没接触过Prepy建议先补前面四篇直接从这篇切入会有点吃力。1. CAMS-Inversion数据到底解决什么问题1.1 区域模拟的浓度边界不是小事区域大气化学模拟有个基本事实模拟域是“开口”的侧边界每小时都在进空气。模拟域内的CO₂浓度变化不只是本地排放和植被呼吸的贡献更大程度上是被“从外面进来的空气”稀释或抬升的。如果边界处给的CO₂浓度是错的哪怕本地通量算得再准出来的浓度场也会整体偏移而且偏移量会随风场一路传播到下风向站点。这就解释了为什么做WRF-GHG模拟时边界条件不能拍脑袋写成固定值。WRF-GHG的碳浓度初始化和侧边界需要一套能反映真实大尺度浓度梯度的全球数据来做驱动。CAMS-Inversion就是目前用得最顺手的全球反演产品之一。1.2 CAMS-Inversion与普通排放清单数据的区别很多人第一次接触CAMS时容易把它的排放清单场和反演浓度场搞混。CAMS本身是数据集系列既提供“自下而上”的排放清单也提供“自上而下”的同化反演结果。边界条件要的是后者也就是经过数据同化系统融合了卫星、地面站点观测后的全球浓度三维场。这里的关键在于“反演”这两个字。反演产品的浓度场不是纯模式外推而是受到观测约束的所以它在全球背景浓度、季节趋势和区域尺度的梯度上都比单纯的气候态或逐月平均更接近真实大气。WRF-GHG-Prepy工具里处理CAMS-Inversion数据目的就是把这种反演浓度场插值到WRF网格上生成wrfbdy侧边界和wrfinput初始浓度。1.3 作为边界条件它与气象场的主从关系不少初学者会问一个问题CAMS数据本身不是也带气象场吗能不能直接用它的气象场替代GFS/FNL理论上同名同源但实际上不建议这么干。CAMS反演的同化系统更关注浓度演变其气象场的内在一致性虽然没问题但分辨率、更新时次和WRF初始化的匹配程度都不如专用气象驱动场。Prepy工具的正常路径是气象场走WPS的met_em温室气体浓度场走CAMS-Inversion或其他全球浓度产品最后在生成IC/BC的时候把两者合并到一起。记住这个主从关系后面看函数代码就不会绕晕。2. WRF-GHG需要什么样的IC/BC文件2.1 wrfbdy与wrfinput里的温室气体变量组织方式WRF-GHG对温室气体的处理本质上是在标准WRF基础上增加了几个额外的三维浓度变量。以CO₂为例real.exe输出的wrfinput里通常有CO₂初始浓度三维场而wrfbdy里则是对应侧边界随时间变化的序列。讲变量命名之前先说一个容易让人崩溃的细节wrfbdy里并不像wrfout那样直接存一个完整的四维CO₂场而是把边界区的点按“四条边”拆分组织。常见的变量形式类似CO2_BT、CO2_BTE这类带侧边标记的名字表示边界上的CO₂浓度随时间的变化。具体命名在WRF-GHG的不同补丁版本里可能有差异但只要理解了“侧边拆分时间序列”的组织方式后面读Prepy的写出函数就顺理成章。2.2 CAMS原始网格与WRF网格的三大差异把CAMS数据转成WRF的IC/BC本质上是在做网格转换难点集中在三个差异上。第一个是水平投影不同。CAMS是规则经纬度网格分辨率在0.5度到1.25度之间WRF是兰伯特正形投影或极射赤面投影的嵌套网格。这个差异需要用水平插值来处理典型的做法是双线性插值如果区域跨度大则要考虑质量守恒。第二个是垂直坐标不同。CAMS是等压面坐标通常在几十层气压层上输出WRF是地形跟随坐标模式层在地面附近加密。两者之间需要把气压面浓度转换到WRF模式的eta层上一般先算出WRF每个格点每一层的实际气压再在垂直方向做对数压力插值。第三个是时间分辨率不同。CAMS反演产品常见是3小时或6小时输出一次而WRF的侧边界需要与模拟时次对齐通常是每小时或每3小时一个时次。中间缺失的时间层靠线性时间插值补齐这个步骤最容易被忽略因为缺失一两层短期看不出问题跑三天以后边界输入的突变会直接影响模拟稳定性。2.3 met_em是中间载体不是替代品Prepy工具生成温室气体IC/BC时通常会读取WPS输出的met_em文件然后把CAMS的浓度场插值到同样的水平和垂直网格上最终写出wrfbdy和wrfinput。很多教程会把这一步讲成“在met_em上叠加浓度”这个描述大方向对但容易让人误解为met_em本身会被修改。实际上met_em文件是不动的它是“模板”和“骨架”。Prepy从中读取网格坐标、气压、地形等几何信息CAMS浓度场经过插值后填入骨架生成新的边界数据。保持met_em原封不动还有个好处real.exe或WRF-GHG启动时可以继续沿用原来的命名和时次逻辑不至于把气象部分弄乱。3. Prepy工具CAMS模块函数解析3.1 模块层级与设计思路注意Prepy具体版本的函数名会迭代不同分支可能有不同的类名。但从源码阅读的常见实现来看CAMS处理链路大体拆成了四段读文件、水平插值、垂直插值、组合写出。我在下面按同一套逻辑结构描述你对照自己的代码找对应函数即可原理是通用的。通常这个模块会定义一两个核心类比如CamscfData用来装载CAMS反演文件的内容。构造函数接收文件路径内部成员变量包括经纬度数组、气压层数组、时间维度和浓度数据。类的核心方法大致分为三类数据读取与变量映射、插值变换、与met_em结合写IC/BC。3.2 数据读取与变量映射read_ncfileCAMS的NetCDF文件直接用netCDF4或xarray都可以打开关键是变量名的映射。原始文件里常见的CO₂浓度变量名可能是co2单位是干空气摩尔分数数值量级通常在390~420 ppm附近。但某些版本会使用质量混合比kg/kg这时就要做单位换算把质量混合比乘以干空气分子量与CO₂分子量的比值再换算成ppm。这一步出错边界浓度能偏差两个数量级整个模拟就废了。读取函数内部还有个容易踩坑的是时间坐标解码。NetCDF里的时间通常是hours since 1900-01-01或seconds since 1970-01-01这类形式如果直接当作数字来处理后面时间插值必然错位。成熟的实现会先统一转成UTC datetime对象再在内存中生成连续时间序列保证与WRF的模拟时次对齐。伪代码大致长这样def read_ncfile(self, ncfile): ds netCDF4.Dataset(ncfile) self.lons ds.variables[longitude][:] self.lats ds.variables[latitude][:] self.pressure ds.variables[lev][:] # 等压面, hPa time_units ds.variables[time].units self.datetimes nc_time_to_datetime(ds.variables[time][:], time_units) co2_raw ds.variables[co2][:] # 如果单位是 kg/kg则换算为 ppmv self.co2 unit_conversion(co2_raw, ds.variables[co2].units) ds.close()3.3 网格转换与插值interp_系列水平插值函数通常接收目标经纬度数组来自met_em的网格和CAMS浓度场输出插值后的三维浓度。实现上多数人会选择scipy.interpolate.RegularGridInterpolator或xarray的interp方法做双线性插值。如果你处理的区域面积很大、跨越多个时区建议先检查CAMS网格与WRF网格的范围避免在数据覆盖区外插值。垂直插值是个更讲究的环节。CAMS的等压面层次从地面一直延伸到高层但WRF模式顶层对应的气压往往远低于CAMS的最低气压层。大多数实现的策略是只做内插不往后外推。也就是说当WRF某一格点的高层气压小于CAMS最顶层气压时直接用CAMS最顶层的浓度值填充而不是强行延伸一段不可靠的曲线。用气压对数坐标做垂直插值通常比线性坐标更平滑尤其在边界层内。注意在低层CAMS的地面气压和WRF的地形气压会有差值如果CAMS数据有地表气压变量最好先读出来作为最低参考层气压避免将浓度插值到地表以下的无效区域。经验不足时最容易在这一步得到一堆NaN后文会专门说这个问题。结合实现垂直插值这一步可以拆成两层循环外层遍历水平格点内层对每个格点做一维垂直插值。遇上大网格和长时间序列这种逐点循环会很慢优化方式是用numpy的apply_along_axis或者直接复用interpolate的轴参数但为写清楚逻辑初版逐点循环更直观。另外还有时间插值函数。它根据当前输出时次找到CAMS前后两个时间层对浓度场做线性时间插值。如果CAMS是3小时间隔而WRF边界是每小时输出时间插值函数就会在每个小时内生成两个中间场。这个函数逻辑简单但执行顺序很有讲究建议先做水平插值再做时间插值最后做垂直插值。原因在于水平插值计算量最大先做时间插值会让参与水平插值的层数翻几倍影响效率。3.4 与met_em结合并写出IC/BCwrite_icbcwrite_icbc是收口的函数。它读取met_em里的目标网格纬度、经度、气压、地形等调用前面的插值函数得到浓度三维场然后写入新的NetCDF文件。初始场部分比较简单直接把第一时次的浓度写入wrfinput对应变量即可。侧边界部分则需要按WRF的边界存储规则做拆分。常见的做法是提取模拟域的四条边每条边各有上下左右不同的索引范围把浓度数据按边写入*_BT、*_BTE等变量。这里最容易错的是边界宽度和索引方向尤其是嵌套域的内边界与外域的重叠区域边界点数量会随bdy_width变化而不同。稳妥起见先摸清你装的WRF-GHG补丁生成了哪些变量名直接参考real.exe输出模板的维度定义来写比想当然更省时间。3.5 关键参数与质量控制从函数解析的角度讲有几个参数对结果影响很大值得花一分钟确认。第一是层顶扩展。如果模拟域高层需要与外部模式衔接最好在垂直插值前先判断CAMS顶层对应的高度必要时对高层浓度做随高度渐变的平滑避免顶边界处出现浓度突变。第二是缺测值处理。CAMS某些版本在极地区域或边缘经度带会有mask如果不填补插值后会出现条带状空洞。简单有效的办法是先用最近邻插值填补缺测再做双线性插值。第三是精度控制。写入NetCDF时用float32比float64更省空间也能与WRF默认精度一致但注意做差值和统计时还是在float64下完成最后写文件前再cast。对新手来说把质量控制流程拆到函数调用前做比在插值后补救要省力。拿到CAMS文件的第一件事应该是画一张全球浓度图看一眼分布是否平滑、量级是否符合预期。这一眼能帮你过滤掉一半以上的数据问题。4. 实操流程从CAMS文件到WPS与REAL4.1 下载与预处理实操第一步是下载对应模拟时段的CAMS-Inversion浓度数据。下载时注意产品类型字段务必确认自己下载的是全球浓度反演场而不是排放清单。时间范围要覆盖整个模拟时间段建议额外多下载前后各一天的数据给时间插值预留缓冲。想起一个常见的坑CAMS数据的发布有时会延迟几天如果你刚跑完今天就要做实时模拟最后几天可能没有最新的反演产品。临时方案是把最后可用时次的浓度场复制填充到后续时次跑业务够用正式研究则建议等数据全了再启动。4.2 运行脚本与配置示例以常见的Prepy调用方式为例运行命令通常由日期、模拟域编号、输入输出路径组成大致长这样python run_prepy.py \ --start 2019-08-01_00:00:00 \ --end 2019-08-05_00:00:00 \ --interval 3 \ --cams-dir /path/to/cams/ncfiles \ --met-dir /path/to/met_em \ --out-dir /path/to/output脚本执行过程中会打印日志有助于定位问题。正常的日志顺序一般是加载CAMS文件列表、读取met_em网格信息、逐时间层插值、合并写文件。如果你看到某一步卡了很久多半是垂直插值循环量比较大属正常现象但如果某一步立刻报错且指向维度不匹配优先检查CAMS和met_em的水平网格范围与垂直层数。4.3 输出检验要点生成的文件不要直接拿去跑real.exe先做两个快速检验。第一个检验是数值范围。用ncdump或Python打开输出文件确认CO₂浓度在300~500 ppm这个合理区间。如果出现上万或负值基本可以断定是单位换算或插值炸了这时回去检查CAMS原始变量的units属性。第二个检验是空间分布。画一张初始时刻的浓度水平分布图重点看边界附近有没有条带状的锯齿纹路。如果边界处浓度和区域内部差异过大说明侧边界变量在写入时可能把索引整体平移了一个格点这个错误跑模式时才会暴露检查起来比较费劲。5. 常见问题与排查技巧实录5.1 单位与数量级错误CAMS数据单位变来变去这是最没有悬念的坑。有的版本用10^-9即ppb级有的用10^-6即ppm级还有的用kg/kg。我处理的早期版本里CO₂浓度看起来在0.0004左右就是质量混合比如果直接当ppm写入结果整个模拟域CO₂浓度等于没有。排查办法很简单打印原始变量的units和数值范围先把单位统一成ppm后再进入插值流程。5.2 时间覆盖缺口与偏移模拟窗口比CAMS数据可用范围长是典型情况。我的处理方式是先检查数据文件的时间列表确认覆盖了模拟窗口。若末尾缺少时次先复制最后时次数据做“保持值”填充跑通流程正式实验要补齐数据。另一个隐蔽问题是时间偏移。CAMS部分产品的时间标签是0点、6点、12点、18点但有些版本使用的是“平均时段中心”或“预报初始时刻”整体偏移3小时或6小时。如果边界条件整体偏移一个时次单看浓度分布看不出异样一旦对比站点观测会发觉模拟浓度的日变化相位整体错位。排查时用站点观测或背景浓度曲线做个简单比对即可。5.3 插值后出现NaN或极端值NaN通常出现在垂直插值阶段。WRF地形很高的格点其地面气压远低于海平面气压如果CAMS数据集没有地表气压变量插值器可能试图在低于最低参考面处找值自然得到NaN。解决路径是读取CAMS的地表气压或地面高度把WRF格点气压低于可插值范围的所有层都统一填成该格点最底层浓度即最近邻向下填充。极端值往往来自水平插值时的mask边界。CAMS在经纬度边缘有大量NaN直接用双线性插值会把NaN扩散到附近格点。建议先对原始浓度场做一层边界填充再进入插值效果立竿见影。5.4 与VPRM通量时间一致性最后一个易错点跟VPRM有关。CAMS反演浓度场的时间体系是UTC而VPRM通量驱动数据用的是MODIS过境时刻和本地时间处理不好会在日尺度上产生时区错位。实际操作中生成边界条件的时间轴统一用UTC不要混入本地时间VPRM通量在写入WRF时已有独立的时间插值逻辑两者分开处理互不干扰。5.5 嵌套域边界不连贯三层嵌套模拟时内层域的外边界要从外层域结果里插值得到。如果直接用CAMS数据分别生成两个域的独立边界可能出现内外层浓度在衔接处跳变。我在实际项目中习惯先只生成最外层域的IC/BC跑出wrfout后再用ndown或real的嵌套逻辑生成内层域边界而不是对每个域都独立跑CAMS插值。这样既省时间也避免了边界不连续带来的数值噪声。最后说点实在的经验CAMS-Inversion处理这个环节我踩过的坑比写代码的时间还多。最初几次跑出来的浓度场整体系统性偏低排查了两三天才发现是单位换算里少乘了分子量比。后来养成一个习惯拿到任何版本的CAMS文件第一件事不是急着跑Prepy而是用xarray打开文件把全球浓度分布和量级先画出来看一眼。这一眼能节省后续一整天。再来一个个人觉得特别重要的小技巧生成IC/BC后先做一次24小时短模拟再对比一下模拟域外围站点或背景站点的CO₂浓度。边界条件的影响在模拟的头一两天最明显如果短模拟的浓度背景和观测趋势一致基本可以放心跑长周期如果不一致优先回去查时间偏移和垂直插值而不是折腾VPRM参数。这个顺序能少走很多弯路。这套处理流程我目前一直在用配合自建的WRF-GHG-Prepy流程从CAMS数据下载到最终wrfbdy生成一次跑通不再返工。后续如果有新的CAMS版本更新或者Prepy工具的函数结构发生变化再看情况补充新文章。